43
edits
No edit summary |
No edit summary |
||
Line 55: | Line 55: | ||
<pre> | <pre> | ||
[mesh] = msh2m_gmsh("fiume","scale",1,"clscale",.1); | [mesh] = msh2m_gmsh ("fiume","scale",1,"clscale",.1); | ||
[mesh] = bim2c_mesh_properties(mesh); | [mesh] = bim2c_mesh_properties (mesh); | ||
</pre> | </pre> | ||
Line 68: | Line 68: | ||
[[File:fiume_msh.png]] | [[File:fiume_msh.png]] | ||
<b> Set the coefficients for the problem:</b> | <b> Set the coefficients for the problem:</b> | ||
Line 83: | Line 82: | ||
<pre> | <pre> | ||
nelems = columns(mesh.t); | nelems = columns (mesh.t); | ||
nnodes = columns(mesh.p); | nnodes = columns (mesh.p); | ||
</pre> | </pre> | ||
Line 95: | Line 94: | ||
<pre> | <pre> | ||
AdvDiff = bim2a_advection_diffusion(mesh, epsilon, 1, phi); | AdvDiff = bim2a_advection_diffusion (mesh, epsilon, 1, 1, phi); | ||
Mass = bim2a_reaction(mesh,delta,zeta); | Mass = bim2a_reaction (mesh,delta,zeta); | ||
b = bim2a_rhs(mesh,f,g); | b = bim2a_rhs (mesh,f,g); | ||
A = AdvDiff + Mass; | A = AdvDiff + Mass; | ||
</pre> | </pre> | ||
Line 103: | Line 102: | ||
<b> To Apply Boundary Conditions, partition LHS and RHS</b> | <b> To Apply Boundary Conditions, partition LHS and RHS</b> | ||
The tags of the sides are assigned by gmsh | The tags of the sides are assigned by gmsh we let <math> \Gamma_D </math> be composed by sides 1 and 2 | ||
and <math> \Gamma_D </math> be the rest of the boundary | |||
<pre> | <pre> | ||
GammaD = bim2c_unknowns_on_side(mesh, [1 2]); | GammaD = bim2c_unknowns_on_side (mesh, [1 2]); ## DIRICHLET NODES LIST | ||
GammaN = bim2c_unknowns_on_side(mesh, [3 4]); | GammaN = bim2c_unknowns_on_side (mesh, [3 4]); ## NEUMANN NODES LIST | ||
GammaN = setdiff (GammaN, GammaD); | |||
jn = zeros(length(GammaN),1); | |||
ud = 3*xu; | jn = zeros (length (GammaN),1); ## PRESCRIBED NEUMANN FLUXES | ||
ud = 3*xu; ## DIRICHLET DATUM | |||
Omega = setdiff (1:length (uin), union (GammaD, GammaN)); ## INTERIOR NODES LIST | |||
</pre> | </pre> | ||
<pre> | <pre> | ||
Add = A( | Add = A(GammaD, GammaD); | ||
Adn = A( | Adn = A(GammaD, GammaN); ## shoud be all zeros hopefully!! | ||
Adi = A( | Adi = A(GammaD, Omega); | ||
And = A( | And = A(GammaN, GammaD); ## shoud be all zeros hopefully!! | ||
Ann = A( | Ann = A(GammaN, GammaN); | ||
Ani = A( | Ani = A(GammaN, Omega); | ||
Aid = A(Ilist, | Aid = A(Ilist, GammaD); | ||
Ain = A(Ilist, | Ain = A(Ilist, GammaN); | ||
Aii = A(Ilist, | Aii = A(Ilist, Omega); | ||
bd = b( | bd = b(GammaD); | ||
bn = b( | bn = b(GammaN); | ||
bi = b( | bi = b(Omega); | ||
</pre> | </pre> | ||
Line 140: | Line 137: | ||
<pre> | <pre> | ||
temp = [Ann Ani ; Ain Aii ] \ [ Fn+bn-And*ud ; bi-Aid*ud]; | temp = [Ann Ani ; Ain Aii ] \ [ Fn+bn-And*ud(GammaD) ; bi-Aid*ud(GammaD)]; | ||
u(Nlist) = temp(1:numel (GammaN)); | |||
u(Omega) = temp(length(un)+1:end); | |||
</pre> | </pre> | ||
<b> Compute the fluxes through Dirichlet sides</b><br> | <b> Compute the fluxes through Dirichlet sides</b><br> | ||
<pre> | <pre> | ||
Fd = Add | Fd = [Add Adi Adn] * u([GammaD; Omega; GammaN]) - bd; | ||
</pre> | </pre> | ||
<B> Compute the gradient of the solution </B> | <B> Compute the gradient of the solution </B> | ||
<pre> | <pre> | ||
[gx, gy] = bim2c_pde_gradient(mesh,u); | [gx, gy] = bim2c_pde_gradient (mesh, u); | ||
</pre> | </pre> | ||
<B> Compute the internal Advection-Diffusion flux</B> | <B> Compute the internal Advection-Diffusion flux</B> | ||
<pre> | <pre> | ||
[jxglob,jyglob] = bim2c_global_flux(mesh,u, | [jxglob, jyglob] = bim2c_global_flux (mesh, u, epsilon, 1, phi); | ||
</pre> | </pre> | ||
<B> Save data for later visualization</B><BR> | <B> Save data for later visualization</B><BR> | ||
The resut can be exported to vtk format to visualize with [[http://www.paraview.org|paraview]] | |||
or [[https://wci.llnl.gov/codes/visit/|visit]] | |||
<pre> | <pre> | ||
fpl_vtk_write_field ("vtkdata", mesh, {}, {[gx; gy]', "Gradient"}, 1); | fpl_vtk_write_field ("vtkdata", mesh, {}, {[gx; gy]', "Gradient"}, 1); | ||
</pre> | </pre> |
edits