I am currently trying to model the cooling of a cylinder via radiation as a function of time using the heat equation, where the emissivity of the material
epsilon_eff is dependent on the temperature. However,
epsilon_eff seems not to be updating for each time step. The basic code is below:
from fenics import * from dolfin import * # Define the domain r1 = 1.0 # inner radius r2 = 2.0 # outer radius L = 3.0 # length N = 20 # number of mesh elements # Create the mesh mesh = CylinderMesh(Point(0, -L/2, 0), Point(0, L/2, 0), r1, r2, N) # Define the function space V = FunctionSpace(mesh, 'P', 1) # emissivity function emiss_in = [0.03,0.042,0.054,0.062,0.080,0.091,0.109,0.116,0.135,0.194,0.213,0.232,0.251] temp_emiss = [300,400,500,600,700,800,900,1000,1100,1600,1800,2000,2200] emiss_f = UnivariateSpline(temp_emiss, emiss_in, k = 2, ext = 3) # Define the initial temperature u0 = Expression('1000', degree=1) # Define the boundary conditions bc = DirichletBC(V, u0, 'on_boundary') # Define the thermal conductivity and radiation coefficient k = Constant(1.0) sigma = Constant(5.67e-8) # Define the variational problem u = TrialFunction(V) v = TestFunction(V) f = Constant(0.0) emiss_eff_array = emiss_f(u.vector().get_local()) emiss_space_eff = FunctionSpace(mesh, 'CG',1) epsilon_eff = Function(emiss_space_eff) epsilon_eff.vector()[:] = emiss_eff_array a = k * inner(grad(u), grad(v)) * dx L = f * v * dx - epsilon_eff * sigma * u**4 * v * ds # Define the time-stepping method T = 100.0 # total simulation time dt = 0.1 # time step t = 0.0 # initial time u = Function(V) # Create the file for storing the solution file = File('cooling.pvd') # Time-stepping loop while t < T: print(t) t += dt solve(a == L, u, bc) print(max(u.vector())) print(max(epsilon_eff.vector())) file << (u, t)
I think that the issue maybe with the line,
epsilon_eff.vector()[:] = emiss_eff_array but I am unsure about how this can be fixed.
Thanks in advance