Take the 2-minute tour ×
Mathematics Stack Exchange is a question and answer site for people studying math at any level and professionals in related fields. It's 100% free, no registration required.

Well, my question is rather from applied maths area, not pure mathematics, so I am not sure that it's a place on this board for one.

I want to solve a direct gravity gradiometry problem on 3D triangular mesh in the $[-1,1]^3$ cube.

Let the area $z>0$ be overground, $z<0$ --- underground. Now push a solid object into the ground and see, how changed the gravity field above the ground.

According to the divergence theorem, we have the field equation: $\nabla \bar{G} = \rho,$ where $\rho$ is the density of mass.

Next, we can build mesh. I used gmsh, my .geo file:

Point(1) = {-1, -1, 0, 1.0};
Point(2) = {-1, 1, 0, 1.0};
Line(1) = {2, 1};
Extrude {2, 0, 0} {
  Line{1};
}
Extrude {0, 0, 1} {
  Surface{5};
}
Extrude {0, 0, -1} {
  Surface{5};
}
Point(25) = {-0.2, -0.2, -0.4, 1.0};
Point(26) = {-0.2, 0.2, -0.4, 1.0};
Line(50) = {25, 26};
Extrude {0.4, 0, 0} {
  Line{50};
}
Extrude {0, 0, -0.4} {
  Surface{54};
}

Physical Volume(84) = {1};  // overground
Physical Volume(85) = {2};  // underground
Physical Volume(86) = {3};  // brick

So, after building mesh and choosing the test functions class, we can build a weak formulation of problem:

$$\int_\Omega \nabla \bar{G} \cdot \bar{v} = \int_\Omega \rho \cdot \bar{v}$$

Now we can run a fem solver. I used sfepy, my .py file:

import numpy as nm

filename_mesh = 'my.mesh'

regions = {
    'Omega'       : ('all', {}),
    'Overground'  : ('nodes of group 1', {}),
    'Underground' : ('nodes of group 2', {}),
    'Brick'       : ('nodes of group 3', {}),
}

field_1 = {
    'name' : 'gravity',
    'dtype' : nm.float64,
    'shape' : (3,),
    'region' : 'Omega',
    'approx_order' : 1,
}

variables = {
    'G'       : ('unknown field',   'gravity',  0 ),
    'g'       : ('test field',      'gravity', 'G'),
}

ebcs = {
}

materials = {
    'm' : ({'rho': {
              'Overground':  1.0e-7,
              'Underground': 1.0e+0,
              'Brick':       1.0e+5
           }},
    ),
    'n' : ({'G' : 1.0 }, )
}

equations = {
    'Gravity' : """dw_div_grad.1.Omega( g, G ) = dw_div.1.Omega( m.rho, g )"""
}

solvers = {
    'ls' : ('ls.scipy_direct', {}),
    'newton' : ('nls.newton', {
        'i_max'      : 1,
        'eps_a'      : 1e-10,
    }),
}

The problem is that the results a little bit unadequate. I expected the field, uniformly directed to the brick, but got that:

enter image description here

What I have to do to correct the model?

Thanks in advance.

share|improve this question
add comment

1 Answer 1

up vote 1 down vote accepted

The first problem I see is that the mesh is wrong - the little box inside is disconnected from the big box it is in:

mesh inside

(The mesh was created by running: "gmsh -3 mesh.geo -format mesh -o mesh.mesh".)

share|improve this answer
add comment

Your Answer

 
discard

By posting your answer, you agree to the privacy policy and terms of service.

Not the answer you're looking for? Browse other questions tagged or ask your own question.