### initial variables values in the whole mesh in a reactional system

37

views

0

Hello,

Currently, I'm trying to integrate a system of PDEs which simulate the chemical species concentration at the end of the mesh. The system is composed of 3 chemical species and current is applied in a 1D domain. The data I receive (from another simulation) is a vector of the initial concentration of the chemical species every couple seconds (I have the value of the elapsed time) at the entry point of the mesh. I want to run the system in such way that I will be able to update the concentrations of the chemical species in the whole mesh prior to solving the system and returning the same structure for a posterior update. So far, I have the following code working properly:

However, my goal is to "fill" the whole mesh with the discrete values obtained in the previous simulation.

How could I manage that? How to assign u and u_n properly (with an array)?

Sincerely,

Andre

Currently, I'm trying to integrate a system of PDEs which simulate the chemical species concentration at the end of the mesh. The system is composed of 3 chemical species and current is applied in a 1D domain. The data I receive (from another simulation) is a vector of the initial concentration of the chemical species every couple seconds (I have the value of the elapsed time) at the entry point of the mesh. I want to run the system in such way that I will be able to update the concentrations of the chemical species in the whole mesh prior to solving the system and returning the same structure for a posterior update. So far, I have the following code working properly:

`from `

fenics` import *`

`T = 5.0 # final time`

`num_steps = 5 # number of time steps`

`dt = T / num_steps # time step size`

`K = 10.0 # reaction rate`

`# Read mesh from file`

`mesh = UnitIntervalMesh(5)`

`# Define function space for system of concentrations`

`P1 = FiniteElement('P', interval, 1)`

`element = MixedElement([P1, P1, P1])`

`V = FunctionSpace(mesh, element)`

`# Define test functions`

`v_1, v_2, v_3 = TestFunctions(V)`

`# Define functions for velocity and concentrations`

`u = Function(V)`

`u_n = Function(V)`

`# Split system functions to access components`

`u_1, u_2, u_3 = split(u)`

`u_n1, u_n2, u_n3 = split(u_n)`

`# Define source terms`

`f_1 = Constant(1)`

`f_2 = Constant(1)`

`f_3 = Constant(0)`

`# Define expressions used in variational forms`

`k = Constant(dt)`

`K = Constant(K)`

`w = Constant(1)`

`# Define variational problem`

`F = ((u_1 - u_n1) / k)*v_1*dx + w*u_1.dx(0)*v_1*dx + K*u_1*u_2*v_1*dx \`

`+ ((u_2 - u_n2) / k)*v_2*dx + w*u_2.dx(0)*v_2*dx + K*u_1*u_2*v_2*dx \`

`+ ((u_3 - u_n3) / k)*v_3*dx + w*u_3.dx(0)*v_3*dx - K*u_1*u_2*v_3*dx + K*u_3*v_3*dx \`

`- f_1*v_1*dx - f_2*v_2*dx - f_3*v_3*dx`

`# Create VTK files for visualization output`

`vtk1 = File('reaction_system/u_1.pvd')`

`vtk2 = File('reaction_system/u_2.pvd')`

`vtk3 = File('reaction_system/u_3.pvd')`

`t = 0`

`for n in range(num_steps):`

`# Update current time`

`t += dt`

`# Solve variational problem for time step`

`solve(F == 0, u)`

`# Save solution to file (VTK)`

`_u_1, _u_2, _u_3 = u.split()`

`vtk1 << (_u_1, t)`

`vtk2 << (_u_2, t)`

`vtk3 << (_u_3, t)`

` # Update previous solution`

` u_n.assign(u)`

Under the tutorial, it says that is possible to update the initial conditions at t=0 by replacing u_0 and projecting:`u_0 = Expression(('sin(x[0])', 'cos(x[0]*x[1])', 'exp(x[1])'), degree=1)`

`u_n = project(u_0, V)`

However, my goal is to "fill" the whole mesh with the discrete values obtained in the previous simulation.

How could I manage that? How to assign u and u_n properly (with an array)?

Sincerely,

Andre

Community: FEniCS Project

Please login to add an answer/comment or follow this question.