elastic wave Log Out | Topics | Search
Moderators | Register | Edit Profile

FlexPDE User's Forum » User Postings » elastic wave « Previous Next »

Author Message
Top of pagePrevious messageNext messageBottom of page Link to this message

xu yi (xhfzjlx)
New member
Username: xhfzjlx

Post Number: 2
Registered: 10-2008
Posted on Tuesday, October 28, 2008 - 03:48 am:   

I want to know how to deal with the propagation of acoustic waves in soil(elastic wave).Can anyeone tell me?Thanks very much.
Top of pagePrevious messageNext messageBottom of page Link to this message

xu yi (xhfzjlx)
Junior Member
Username: xhfzjlx

Post Number: 3
Registered: 10-2008
Posted on Thursday, November 20, 2008 - 07:26 pm:   

Can someone show me an example?
Especially how to define the boundary conditions and set the load?
Top of pagePrevious messageNext messageBottom of page Link to this message

xu yi (xhfzjlx)
Member
Username: xhfzjlx

Post Number: 4
Registered: 10-2008
Posted on Thursday, November 20, 2008 - 08:35 pm:   

Problem;
1. How to set the field source? Is it right to load the stress?Is there some other method to set the source load?
2. The top of the region is free surface.I define it with stress and set to zero(value(varialbe)=0 in boundaries),is it right?
3. In the contour plot of x-displacement or Y-displacement,near the top surface,I can't see the rayleigh wave,why?
application/octet-streampropagation of elastic wave in solid
elasticwave.pde (4.1 k)
Top of pagePrevious messageNext messageBottom of page Link to this message

xu yi (xhfzjlx)
Member
Username: xhfzjlx

Post Number: 5
Registered: 10-2008
Posted on Thursday, November 20, 2008 - 08:36 pm:   

Here is my script:

{Analysis the propagation characteristic of semi-infinite soil by FlexPDE.
According the elastic wave theroy, the wave 2D equation is:

(lamda+G)*(dxx(Ux)+dxy(Uy))+G*DIV(GRAD(Ux))=rho*dtt(Ux)
(lamda+G)*(dxy(Ux)+dyy(Uy))+G*DIV(GRAD(Uy))=rho*dtt(Uy)

where,
Ux ---dispalcement of x coordinate(horizon direction)
Uy----dispalcement of y coordinate(vertical direction)
lamda ----Lame's constant,
lamda=nu*E/((1+nu)*(1-2*nu))
nu is the Poisson's Ratio, E is Young's Modulus for Material.
G--shear Modulus for Material
G=E/(2*(1+nu))
rho--density for material

for convenience, introducing two potential function phi and psi satisfied,
Ux=dx(phi)+dy(psi)
Uy=dy(phi)-dx(psi)

The wave equation can be transformed to another form,

dtt(phi)=((lamda+2*G)/rho)*DIV(GRAD(phi))=Vp*Vp**DIV(GRAD(phi));
dtt(psi)=(G/rho)*DIV(GRAD(psi))=Vs*Vs*DIV(GRAD(psi));

The vertical Stess of material is
sigma=2*G*dy(Uy)+lamda*(dx(Ux)+dy(Uy))

The shear Stress is
tao=G*(dx(Uy)+dx(Ux))

The second time derivative in the acceleration term cannot be modelled
directly in FlexPDE, but the problem can still be solved.For example,

VARIABLES u,v
EQUATIONS
U: v = dt(u)
V: c^2*div(grad(u)) = dt(v)


}

TITLE 'wave propagation in soil ' { the problem identification }

select
{ deltat=1.0e-4 { Start out at careful timestep, it will grow. }}
ngrid=30 { Grid a little more densely than default }

regrid = off


variables { Recall that the declared variable range, if too large,
will affect the interpretation of error, and thus the
timestep and solution accuracy }

phi (threshold=1e-7) { Potentials }
psi (threshold=1e-7)
vphi (threshold=1e-5) { Velocities }
vpsi (threshold=1e-5)
sigma(threshold=1e-7)
tao(threshold=1e-7)

definitions
L1=200 {the length of analysis square}
L2=0.06 {force driving zonge}
W1=0.5*L1-0.5*L2
W2=0.5*L1+0.5*L2

Vp= 1800 { Velocity of pulse wave m/s }
Vs= 1260 {Velocity of shear wave m/s }
rho= 1530 {density for soil}

G=Vs*Vs*rho; {here G is mu of lame constants}
lamda=Vp*Vp*rho-2*G;

{dispalcements}
Ux=dx(phi)-dy(psi)
Uy=dy(phi)+dx(psi)
U = sqrt(Ux^2+Uy^2) {displacement amplitude }

{strain}


freq =120 { driving frequency }
period = 1/freq

force = 1 { loading force in Newtons (~1 pound force) }
{ the driving force is sinusoidal in time: }
jiggle =-force*sin(2*pi*freq*t)*upulse(t-0,t-1*period)


initial values
phi=0
psi=0
vphi =0
vpsi =0
sigma=0
tao=0

equations
phi: vphi=dt(phi)
vphi: Vp^2*div(grad(phi))=dt(vphi)
psi: vpsi=dt(psi)
vpsi:Vs^2*div(grad(psi))=dt(vpsi)
sigma: 2*G*dy(Uy)+lamda*(dx(Ux)+dy(Uy))=sigma
tao: G*(dx(Uy)+dy(Ux))=tao

boundaries
region 1
start (0,0)
value (sigma) =0 ! stress boundary condition
value(tao)=0
line to (W1,0)
value(sigma)=jiggle {Apply oscillatory vertical load in }
value(tao)=0
line to (W2,0)
value(sigma)=0
value(tao)=0
line to(L1,0)
natural(phi)=0
natural (psi)=0
line to(L1,-L1)
natural(phi)=0
natural (psi)=0
line to(0,-L1)
natural(phi)=0
natural (psi)=0
line to Close



time 0 to 1000*period

plots
for t= 0 by period/10 to endtime
grid (x,y)
contour(Ux) PAINTED as "X-displacement(M)"
contour(Uy) PAINTED as "Y-displacement(M)"
contour(dt(sqrt(Ux*Ux+Uy*Uy))) PAINTED as "X-Velocity(M/s)"
contour(phi) PAINTED as "P Wave"
contour(psi) PAINTED as "S Wave"
history(dt(Uy)) at (L1/2,-0.1*L1) (L1/2,-0.2*L1) (L1/2,-0.3*L1)
history(dt(Ux)) at (W2,0) (W2+1,0)
{contour(U) PAINTED as "Sound Field"}
contour(sqrt(phi^2+psi^2)) PAINTED as "Sound Field"


end
Top of pagePrevious messageNext messageBottom of page Link to this message

Marek Nelson (mgnelson)
Moderator
Username: mgnelson

Post Number: 84
Registered: 07-2007
Posted on Friday, November 21, 2008 - 12:20 am:   

We are not familiar with the model that you are constructing and so we cannot comment on the meaning of the equations and boundary conditions. These are issues that you will have to determine for yourself.

There are a couple of issues about how to pose your problem to FlexPDE that I can comment on:

1)

It appears as if your acoustic waves are bouncing of the sides of your domain. This is implies that the boundary conditions are not correct.

The NATURAL boundary condition assigns values to the surface integrals generated by integration by parts of second order spatial derivatives (divergence theorem). In your problem, this means :

NATURAL(PHI) and NATURAL(PSI) are undefined.
NATURAL(VPHI) provides the value of (Vp^2)*NORMAL(grad(PHI))
NATURAL(VPSI) provides the value of (Vs^2)*NORMAL(grad(PSI))
NATURAL(SIGMA) provides the value of (2*G+lambda)*NORMAL(Uy) on horizontal surfaces and lambda*NORMAL(Ux) on vertical surfaces
NATURAL(TAO) provides the value of G*Uy on vertical surfaces and G*Ux on horizontal surfaces

See NATURAL in the Help index.

2)

Your thresholds seem high. If I plot your variables, I see that they are in the range of 1e-30. This means that the FlexPDE error control is ignoring most variations in the variables. See THRESHOLD in the Help index.

I tried using a THRESHOLD that is on the same order of magnitude as the values observed, but the time step control cut the time step drastically. This implies a discontinuous initial condition which creates high frequency transients. You should make sure that your initial conditions are self consistent.

Add Your Message Here
Post:
Username: Posting Information:
This is a private posting area. Only registered users and moderators may post messages here.
Password:
Options: Enable HTML code in message
Automatically activate URLs in message
Action:

Topics | Last Day | Last Week | Tree View | Search | Help/Instructions | Program Credits Administration