Just to inform you. Actually I tried to find out what error could lead to non-solving the 3 eqn I had before.
Therefore I just used your navier_stokes.py script with a different mesh, the mesh I also use for my equations and the definitions of the regions and boundary conditions. With the 3d_elbow mesh it worked fine but using the cube_big_tetra.mesh I get the same error like in my script, umfpack out of memory.
Could the reason be somewhere in the mesh/region/boundary condition defintion?
Quinta-feira, 23 de Maio de 2013 15:26:20 UTC+1, Robert Cimrman escreveu:
On 05/23/2013 04:12 PM, dj....@gmx.de javascript: wrote:
Yes, exactly, if we consider incompressiblity in a homogeneous case then both will be uncoupled. That´s what of course works fine to model. But I really want the coupled case. So if I use your proposition, the last equation:
Laplace(V) + (grad(V)/grad(h)) * Laplace (h) = 0
How would I implement it in sfepy? Because this is one equation but 2 unknown variables. I can only use one test variable. Is it possible to solve?
Well, we still have a missing equation. Otherwise the above equation could be implemented for example using a solution-dependent material coefficient, as in [1].
r.
[1]
http://sfepy.org/doc-devel/examples/diffusion/poisson_field_dependent_materi...
Quinta-feira, 23 de Maio de 2013 15:00:53 UTC+1, Robert Cimrman
On 05/23/2013 03:29 PM, dj....@gmx.de javascript: wrote:
Ok, the complete equations are basically those, from which I derived
simplified ones:
- j=-sigma * grad(V) - L * grad(h)
- div j = 0
- u= -k/eta * grad (h)
- L/sigma=grad(V)/grad(h)
With L, sigma, k, eta constant as this is a homogeneous case. Also it is stationary. h is the hydraulic potential (scalar) and V the Self-Potential (scalar) , u is the velocity vector. When I apply eq. 2) on eq. 1) I get the following:
- sigma * Laplace(V) + L * Laplace (h) = 0
This can be modelled with dw_laplace.
Eq. 4) I tried with the dw_stokes but as you are right I thought it takes the gradient, so I will try your proposittion with dw_v_dot_grad_s. If considering a divergence-free velocity field, we could add div u = 0, the incompressibilty,which will lead to a modified eq. 3)
3*) -k/eta * Laplace (h) = 0
Well, then 5) reduces to sigma * Laplace(V) = 0, which certainly can be solved, but I am not sure it is what you want :) - there is no coupling at all.
But I actually didn´t want this case because I tried to add a source (a well) into the middle of the cube. Anyway, for testing purpose we could try to use this. So are eqn. 4), 5), 3*) correct and can they be the 3 eqn. needed to model? Or should I use div u = 0 instead of 3*)?
I know nothing about this problem, so all I can provide are negative statements ("this will not work") - I do not know what is the right weak
3*) means that V and h are solutions to Laplace equations and are decoupled. Then I do not know, if/how 4) could be enforced. Is it that L, sigma actually depend on the solution, so that 5) could be written as this:
Laplace(V) + (L / sigma) * Laplace (h) = 0
so
Laplace(V) + (grad(V)/grad(h)) * Laplace (h) = 0
- is grad(h) always nonzero?
r.
Thanks for your help!
Quinta-feira, 23 de Maio de 2013 12:44:56 UTC+1, Robert Cimrman escreveu:
I think we first need to get the equations right.
- sigma * Laplace(V) + L * Laplace (h) = 0
- L/sigma * grad(h) - grad(V) = 0
- the 2) uses grad(h), grad(V), but you use dw_stokes, which
corresponds
to \int div(h) q, \int \div{V} q (and thus has a scalar test function). Look at dw_v_dot_grad_s instead(?)
Do you have the complete equations written somewhere? Maybe you are over-simplifying something here.
It seems to me, that you have:
- scalar unknown hydraulic head h
- scalar unknown self_Potential signal V
- velocity vector v
So there should be three equations, two of them multiplied by a
scalar
test functions corresponding to h and V, and one multiplied by a vector test function, corresponding to v.
r.
On 05/23/2013 12:48 PM, dj....@gmx.de javascript: wrote:
Thanks Robert, when I try now the corrected equations (scalar test variable q for
escreveu: the formulation. the
first equation, vectortest variable s for the second) and introduce a third equation (simplified Navier Stokes without stokes term, with a vector test variable uv) I get the following memory error:
failed with UMFPACK_ERROR_out_of_memory Both vector test variables are related to the same field (velocity). I might not have understood the formulation of the weak eqn.correctly.
'Head-SP' : """dw_laplace.i3.Omega( fluid.sigma, q, v ) +
dw_laplace.i3.Omega( fluid.C_sigma_rel, q, h ) = 0""", 'Coupling-Coef' : """dw_stokes.i3.Omega( fluid.C_sigma_rel, s, h ) - dw_stokes.i3.Omega( s, v ) = 0""", 'Navier-Stokes' : """dw_div_grad.i3.Omega( fluid.dyn_viscosity, uv , u ) + dw_convect.i3.Omega( uv , u ) = 0""",
I'm using linux fedora as a virtual box in windows 7 system, allowing 4GB of RAM for linux.
Sorry for all this, just need this step and then I hope I will get it. Thanks again, Djamil
Quarta-feira, 22 de Maio de 2013 16:33:25 UTC+1, Robert Cimrman escreveu: > > Hi Djamil, > > On 05/22/2013 04:16 PM, dj....@gmx.de javascript: wrote: >> Dear friends, >> I'm trying to solve the following *equations *for the self-potential >> distribution along a water injection in a well. >> >> 1) sigma * Laplace(V) + L * Laplace (h) = 0 >> 2) L/sigma * grad(h) - grad(V) = 0 > > The weak form of your equations below seems to be wrong - as grad(h) and > grad(V) are vectors, the second equation has to be tested by a vector test > variable. Then this vector test variable should have a corresponding > unknown > vector variable (the velocity) - there should be IMHO one more equation > with that. > > Also, each one of the above equations should be tested by a single test > variable - you are mixing it in the equations. > >> *Here h is the hydraulic head (pressure) distribution, V is the >> Self_Potential signal*. >> I'm using the *3d_cube_big_tetra mesh*. The code is attached. >> >> The *boundary conditions* are: >> Left and right side of the cube h=100m >> Left side of the cube (reference point) V=0mV. >> Center cells of the cube: h=1000m (injection of water) >> >> In general it should be very easy to find a solution with sfepy FE but I >> don´t get any relation between V and h, although both should be coupled > via >> the L/sigma parameter. L and sigma are constant in this *homogeneous > case*(L=-1.5e-7, sigma=1). >> >> Hydraulic head (h) is calculated correctly, but the SP signal (V) shows > no >> connection to h. Changing the parameters L or sigma has absolutely no >> effect on the solution. It depends only on the boundary conditions. >> >> So I tried to implement it with the following *equations*: >> * >> 1 'SP' : >> """dw_laplace.i3.Omega( fluid.sigma, w, v ) + dw_laplace.i3.Omega( >> fluid.C_sigma_rel, q, h ) = 0""", >> 2 'Coupling' : >> """dw_stokes.i3.Omega( fluid.C_sigma_rel, q, h ) - > dw_stokes.i3.Omega( >> w, v ) = 0""", * >> >> A probe image through the middle of the cube is attached. >> >> So my questions: >> >> 1) Am I using the right equations? I also tried it with dw_div_grad > instead >> of dw_laplace (should be the same, no?) but unfortunately the > calculation >> then doesn't converge. > > Use dw_div_grad for vector variables, and dw_laplace for scalar variables. > >> 2) Why is the coupling between both Laplacian eqn. not given? There > should >> be an effect of L and sigma even without defining eqn 2 (Coupling) but > it >> has no effect. This is not a SP problem but more a general issue because >> there should be a coupling between both unless one part is zero. > > I guess this is caused by the one equation that is IMHO missing... > >> That's why I'm afraid that having problems already with a homegeneous > case >> leads to much more problems with inhomogeneous sigma and L > distributions. >> >> Thanks very much for any hints and suggestions. I'm really out of ideas > now. > > I some more observations/suggestions related to the example file: > > - in field_2, use 'shape' : 'vector' instead of 'shape' : (2,), - that > will > work for both 2D and 3D meshes > - it does not run with current sfepy (fmf_sumLevelsMulF(): ERR_BadMatch (4 > == > 12, 1 == 1, 1 == 1)) - this is caused by using the scalar test functions > in the > second equation > - when I correct the test field of the second equation to be "s" (vector), > the > system assembles, but the matrix is singular, as there is no "u" in the > equations.