Hi Alec,
On 08/21/2012 10:50 AM, Alec Kalinin wrote:
Dear SfePy users,
I am new to the SfePy and first of all I want to say thank you to the developers of this FEM engine. Now I am trying to use SfePy for my tasks and I got several questions.
Consider problem in the example ``examples/diffusion/poisson.py'' (I will use LaTeX notations further for the math). We need to find unknown function $t(x)$ such that $$ c \Delta t = 0, x \in \Omega, $$ $$ t(x) = T_1, x \in \Gamma_{D1}, $$ $$ t(x) = T_2, x \in \Gamma_{D2}, $$ $$ \frac{\partial t(x)}{\partial n} = 0, x \in \Gamma_{N}, $$ where $T_1$ and $T_2$ is known functions, $\Gamma_{D1} = {x \mid x < 0.00001}$, $\Gamma_{D2} = {x \mid x > 0.099999}$, $\Gamma_{N} = \partial \Omega \setminus (\Gamma_{D1} \cup \Gamma_{D2})$. This is the mixed boundary value problem with both Dirichlet and Neumann boundary conditions.
Could you, please, answer me on following questions:
- Is it possible to get a flux on the $\Gamma_{D1}$ and $\Gamma_{D2}$
surfaces after the problem has been solved?
Not yet - we actually never saved a solution on a surface. I will try to fix this issue now. Then, something like the following code should work (see [1] for definition of d_surface_flux term, K is a permeability matrix):
- Define a postprocessing function.
def post_process(out, pb, state, extend=False): """ Calculate fluxes. """ from sfepy.base.base import Struct from sfepy.fem import extend_cell_data
ev = pb.evaluate
flux = ev('d_surface_flux.2.Gamma(m.K, t)', mode='el_avg')
flux = extend_cell_data(flux, pb.domain, 'Gamma') # Extend to whole domain
- this does not work now. out['flux'] = Struct(name='output_data', mode='cell', data=flux, dofs=None) return out
- Add the function to the options:
options = { 'post_process_hook' : 'post_process', }
- Suppose we solve pure Dirichlet problem (no Neumann conditions). We need
to specify essential boundary conditions on the all surface nodes? What is the best way to select all surface nodes?
You can use the following region definition:
regions = { 'Gamma' : ('nodes of surface', {}), }
to select all the surface nodes.
- Suppose the Neumann boundary conditions is not zero. In this case we
need to add surface integral with the flux. Do we need to add new term for the flux in the problem description? And this term will be the surface integral but in the problem description we have only tetrahedral domain 3D mesh. What is the best way to select only surface part of the mesh for the surface integration?
Yes, you need to add the flux term "dw_surface_ndot.2.Gamma(m.c, s)" to the equation, where m.c. is the given flux (a material), and s is the test function, again check [1], c is the flux vector. As for the surface selection - see above.
See also example [2] - linear elasticity with surface tractions.
Cheers, r.
[1] http://sfepy.org/doc-devel/terms_overview.html [2] http://sfepy.org/doc-devel/examples/linear_elasticity/linear_elastic_tractio...