Hi Alexander
On 06/19/2013 04:42 AM, Alexander Mihailov wrote:
Hi all.
I have already asked my question at http://math.stackexchange.com/questions/422928/modeling-gravity-field-with-f... and http://stackoverflow.com/questions/17160835/modeling-gravity-field-with-fini... but seems there are no people there, who knows sfepy and finite elements, in particular.
I want to model the gravity field on 3D rectangular area. It can be described by the equation: div(G) = rho. Here G is vector unknown function, rho is scalar parameter, which is constant at the selected point and fully determined by material.
The weak form is: int( div(G) * g ) = int( rho * g ). Here g is vector test function. Using sfepy syntax: {'Gravity' : """dw_div_grad.1.Omega( g, G ) = dw_volume_lvf.1.Omega( m.rho, g )"""}
The terms you use do not correspond to the equation. As div(G) is a scalar, it cannot be used with a vector test function, but a scalar test function - the term that does that is dw_stokes, but then the matrix would not be square.
Try using a different approach - use the scalar gravity potential u as your function, such that G = grad(u). As div(grad(u)) = rho is the Poisson equation, use dw_laplace and dw_volume_lvf terms with scalar unknown and test fields. Check the equation in [1]. Then G could be computed in a postprocess hook function by evaluating ev_grad term - see [2] for an example of a similar hook function (post_process()).
r.
[1] http://sfepy.org/doc-devel/examples/diffusion/poisson_parametric_study.html [2] http://sfepy.org/doc-devel/examples/biot/biot_npbc.html
I used gmsh to build .mesh file. My .geo, .mesh and .py files attached to this message.
When I try to run sfepy, it says that matrix is singular, so the solving process does not converge.
What I have to do to correct the model?
Thanks in advance.