Gaussian spatial push on 2D surface in 3D geometry Answered
Hi there,
I'm looking for an efficient way to apply mechanical load with gaussian spatial amplitude to avoid edge artifacts for a point source push applied on a 3D anisotropic material (e.g. tendon). I've done this presently by pushing using 5 concentric point sources (cylindrical push) having 5 radii (e.g., 10micron to 50 micron with 10 micron interval). A DO loop was used to apply each push; at each execution the radius of the push cyln enlarges and the amplitude picks a gaussian amplitude reduction factor from a separate file. Since its a surface push applied through DO loop, my method is so slow so that I can't use more resolved amplitude factors to make a smooth push (e.g., by using 50 or 100 radii or circles to mimic a better gaussian scale factor). Anyone has an idea to resolve this?
Thank you
Maju
9 comments
Section of code to make the push..
c Import time data of Gaussian
data
hist gauss_data * supergaussian.dat $push_duration 1.0 $push_offset
end
c Import spatial data of push
data
hist spatial_data * testGaus.dat $push_width_scale $push_amplitude
end
c Push definition
func hist gauss_data 1.0
symb sptsize = 2e-4
symb nsptloop = 5
do push_loop ii 1 $nsptloop 1 /*$i$idx-1 1
symb #get { rxcrd rycrd rzcrd } crdnode $ii 1
symb #get { push_scale } histval spatial_data $rxcrd
plod
symb radc = $ii * $sptsize / $nsptloop /* calculate radius of cyln
chek off off
pdef push_loading$ii func $push_scale
vctr vct$ii 0. 1. 0.
cyln push_loading$ii vct$ii stnd y 0 $grid_spacing 0.0 0.0 $radc
end
end$ push_loop
Hi Maju,
When you say the more resolved simulations are so slow, how slow do you mean? There are a few ways models can be sped up but I would need to take a look at the full code to make suggestions. Are you able to share the files via your preferred cloud sharing service?
Best Regards,
Chloe
To give little more explanation to the above mentioned problems:
Here below is a 3D simulation with grid size 5 micron (earlier mentioned was 15 micron)
Displacement fields after 333 micro seconds for an isotropic (mu = G = 10k) system is displayed in the below figs those spatially filtered at different bandwidths (previously BW = 500 to 10kHz)
Filtering values are in the title : BW [10 10000] = 10Hz to 10kHz.
As you can see, at lower bandwidth, a square pattern is visible in the diagonal directions. At higher cut off of the low limit, squares were almost invisible but the shape (it should be circle) was lost.
Thank you
Maju
code used:
Hi Maju,
We typically suggest creating point sources by adding a circular primitive with a dummy material and loading the outside surface of that primitive using the material interface. So, for your example you could copy the dermis material e.g. dermis2 which will have the same material properties as dermis. In the site command you would add a point source at the desired location using cyln or sphr (depending if 2D or 3D) and give it material dermis2. You can then load the interface between the materials dermis and dermis2 using the command plod sdf2 which would give you your point source. See the command reference for help on the code syntax.
I hope this helps.
Best Regards,
Chloe
Thank you for the prompt response. I'll check it out.
Maju
Hi Chloe,
I've implemented sphere push that resolved the problem of 'slowing down' of my program. Longitudinal wave propagation seemed correct, which shows a nice circular pattern radiating away from the point like source on a 3D geometry with Z axis as depth.
Now, if I wait for shear wave,
then the circular shape of the velocity field (in Z for a push on XY plane in a 3D model) is not preserved. please see the picture taken at 4.24ms. Shear velocity used was 3m/s
c Push
symb push_duration = 50E-6
symb push_amplitude = 5000.0
symb push_offset = 0E-6 /* $push_duration / 2.0 * 5.0
symb sptsize = 50e-6
symb push_width = $sptsize
symb push_width_scale = $push_width
c Materials - Layer
symb density = 1000.0
symb dermis_mu = 10000.0
symb dermis_G = 1 * $dermis_mu
symb q2 = 0
symb dbm = 0 /* delta by mu
symb q1 = 0.5 * ( $q2 - $dbm * $dermis_mu )
symb param_lame_lambda = ( 35.0 * 35.0 - 2.0 ) * $dermis_mu /* $density * $cp * $cp - 2 * $shearmod_subq_mu /*( 1500.0 * 1500.0 - 2.0 ) * $shearmod_subq_mu
symb c11 = $param_lame_lambda + 2.0 * $dermis_mu + $q2
symb c12 = $param_lame_lambda + $q1
symb c22 = $param_lame_lambda + 2.0 * $dermis_mu
symb c23 = $param_lame_lambda
symb cp = sqrt ( $c11 / $density )
symb c44 = $dermis_mu
symb c55 = $dermis_G
symb c66 = $dermis_G
symb dimxy = 0.015
symb dimz = 0.001
symb #keycord x 1 0 $dimxy $dimxy
symb #get { idx } rootmax x
symb #keycord y 1 0 $dimxy $dimxy
symb #get { jdx } rootmax y
symb #keycord z 1 0 $dimz
symb #get { kdx } rootmax z
symb box = 0.00005
symb #keyindx i 1 $idx 1 $box 1
symb #keyindx j 1 $jdx 1 $box 1
symb #keyindx k 1 $kdx 1 $box 1
grid $i$idx $j$jdx $k$kdx
symb domainxy = $dimxy + $dimxy
geom
xcrd 0 $domainxy
ycrd 0 $domainxy
zcrd 0 $dimz
matr
wvsp off
type lean
prop dermis $density * * 0.005
stif $c11 $c12 $c12 0 0 0 $c22
$c23 0 0 0 $c22 0 0
0 $c44 0 0 $c55 0 $c66
c c hrgl subq_layer_matr BB
end
c matr
c wvsp on
c type elas
c prop steel 1000. 110 3
c end
matr
wvsp on
type elas
prop steel2 $density $cp 3
end
site
regn dermis
sphr steel2 $dimxy $dimxy 0 0.0003
end
c Boundary conditions
boun
side xmin absr
side xmax absr
side ymin absr
side ymax absr
side zmin free
side zmax absr
end
grph
plot matr
end
term
data
hist gauss_data * supergaussian.dat $push_duration 1.0 $push_offset
end
func hist gauss_data 1.0
plod
pdef pld1 func $push_amplitude
spot spot1 $dimxy $dimxy 0
sdf2 pld1 spot1 dermis steel2
end
calc
pres acoustic
max aprs none apmx
disp
end
prcs
proc plot save
exec 1
grph
plot aprs
end
end$ proc
proc plot 5
term
symb simtime = 0.005
c
c Runtime calculations
c
symb #get { step } timestep
symb nexec = $simtime / $step
symb nloops = 100
symb nexec2 = $nexec / $nloops
symb #get { tbeg } wtime
do loopi ii 1 $nloops 1 /* setting up do loop
exec $nexec2 /* run model for 'nexec2' timesteps
c
c FLXDATO OUTPUT
c
data
file out skin_multilayer_results_surface$ii
out1 tmas 1 1 1 1
out1 xcrd
out1 ycrd
out1 zcrd
out1 matd
out1 mats
out1 matp
c ---
c To export only surface data, append * * 1 1 to these variables
out1 xvel * * * * 1 1
out1 yvel * * * * 1 1
out1 zvel * * * * 1 1
out1 xdsp * * * * 1 1
out1 ydsp * * * * 1 1
out1 zdsp * * * * 1 1
c out1 xcrl * * * * * *
c out1 ycrl * * * * * *
c out1 zcrl * * * * * *
end
end$ loopi
symb #get { tend } wtime
symb trun = $tend - $tbeg
symb #get { label } jobname
symb #save '$label.symb'
stop
Hi Maju,
I think the absorbing boundaries may be introducing artefacts at the corners. To combat this you can bend your model into a circle so the absorbing boundary is a smooth surface. This will also ensure the point source surface is smooth (instead of stepped due to the mesh). Here is an example that shows you the syntax for the bend command.
Best Regards,
Chloe
Please sign in to leave a comment.
Didn't find what you were looking for?
New post