[cig-commits] [commit] devel,master: fixed the calculation of pressure in the case of an anisotropic material (it is then necessary to compute sigma_yy using additional c_IJ anisotropic parameters); updated the manual and the Par_file accordingly; added a new example: EXAMPLES/anisotropic_zinc_crystal; also removed useless white spaces (ca0ebca)
cig_noreply at geodynamics.org
cig_noreply at geodynamics.org
Wed Jun 18 15:23:41 PDT 2014
Repository : https://github.com/geodynamics/specfem2d
On branches: devel,master
Link : https://github.com/geodynamics/specfem2d/compare/fc67e6fd7ad890705b2b72b4b3c509accb22249e...e9ca46c40131588d89d7b0883250bc6584ce6b4c
>---------------------------------------------------------------
commit ca0ebcafb9ed6b1df46bf70165bd716268215823
Author: Dimitri Komatitsch <komatitsch at lma.cnrs-mrs.fr>
Date: Mon Sep 16 23:48:28 2013 +0000
fixed the calculation of pressure in the case of an anisotropic material (it is then necessary to compute sigma_yy using additional c_IJ anisotropic parameters);
updated the manual and the Par_file accordingly;
added a new example: EXAMPLES/anisotropic_zinc_crystal;
also removed useless white spaces
>---------------------------------------------------------------
ca0ebcafb9ed6b1df46bf70165bd716268215823
.../Par_file | 73 ++++++++++++----------
.../SOURCE | 6 +-
.../topoaniso.dat | 10 +--
3 files changed, 47 insertions(+), 42 deletions(-)
diff --git a/LuoYang_fluid_solid_kernel/Par_file b/anisotropic_zinc_crystal/Par_file
similarity index 77%
copy from LuoYang_fluid_solid_kernel/Par_file
copy to anisotropic_zinc_crystal/Par_file
index 0a2f4ae..db9e10f 100644
--- a/LuoYang_fluid_solid_kernel/Par_file
+++ b/anisotropic_zinc_crystal/Par_file
@@ -1,17 +1,18 @@
+
# title of job
-title = Test
+title = Test of anisotropic zinc crystal
# forward or adjoint simulation
-SIMULATION_TYPE = 1 # 1 = forward, 2 = adjoint + kernels
-NOISE_TOMOGRAPHY = 0 # 0 = earthquake simulation, 1/2/3 = noise simulation
-SAVE_FORWARD = .true. # save the last frame, needed for adjoint simulation
+SIMULATION_TYPE = 1 # 1 = forward, 3 = adjoint + kernels; 2 is purposely UNUSED (for compatibility with the numbering of our 3D codes)
+NOISE_TOMOGRAPHY = 0 # 0 = earthquake simulation, 1/2/3 = noise simulation
+SAVE_FORWARD = .false. # save the last frame, needed for adjoint simulation
# parameters concerning partitioning
nproc = 1 # number of processes
partitioning_method = 3 # SCOTCH = 3, ascending order (very bad idea) = 1
PERFORM_CUTHILL_MCKEE = .false. # perform inverse Cuthill-McKee (1969) optimization/permutation for mesh numbering (can be very costly and not very useful)
-ngnod = 9 # number of control nodes per element (4 or 9)
+ngnod = 4 # number of control nodes per element (4 or 9)
initialfield = .false. # use a plane wave as source or not
add_Bielak_conditions = .false. # add Bielak conditions or not if initial plane wave
assign_external_model = .false. # define external earth model or not
@@ -23,8 +24,8 @@ freq0 = 10 # frequency for viscous attenua
p_sv = .true. # set the type of calculation (P-SV or SH/membrane waves)
# time step parameters
-nt = 1600 # total number of time steps
-deltat = 1.1d-3 # duration of a time step (see section "How to choose the time step" of the manual for how to do this)
+nt = 2000 # total number of time steps
+deltat = 55.e-9 # duration of a time step (see section "How to choose the time step" of the manual for how to do this)
USER_T0 = 0.0d0 # use this t0 as earliest starting time rather than the automatically calculated one
time_stepping_scheme = 1 # 1 = Newmark (2nd order), 2 = LDDRK4-6 (4th-order 6-stage low storage Runge-Kutta), 3 = classical RK4 4th-order 4-stage Runge-Kutta
@@ -37,10 +38,10 @@ N_SLS = 2 # number of standard li
f0_attenuation = 5.196152422706633 # (Hz) relevant only if source is a Dirac or a Heaviside, otherwise it is f0 the dominant frequency of the source in the CMTSOLUTION file
# receiver set parameters for seismograms
-seismotype = 6 # record 1=displ 2=veloc 3=accel 4=pressure 5=curl of displ 6=the fluid potential
+seismotype = 1 # record 1=displ 2=veloc 3=accel 4=pressure 5=curl of displ 6=the fluid potential
NSTEP_BETWEEN_OUTPUT_SEISMOS = 5000000 # every how many time steps we save the seismograms (costly, do not use a very small value; if you use a very large value that is larger than the total number of time steps of the run, the seismograms will automatically be saved once at the end of the run anyway)
save_ASCII_seismograms = .true. # save seismograms in ASCII format or not
-save_binary_seismograms_single = .true. # save seismograms in single precision binary format or not (can be used jointly with ASCII above to save both)
+save_binary_seismograms_single = .false. # save seismograms in single precision binary format or not (can be used jointly with ASCII above to save both)
save_binary_seismograms_double = .false. # save seismograms in double precision binary format or not (can be used jointly with both flags above to save all)
SU_FORMAT = .false. # output single precision binary seismograms in Seismic Unix format (adjoint traces will be read in the same format)
subsamp_seismos = 1 # subsampling of the seismograms to create smaller files (but less accurately sampled in time)
@@ -51,11 +52,11 @@ rec_normal_to_surface = .false. # base anglerec normal to surfa
save_ASCII_kernels = .true. # save sensitivity kernels in ASCII format (much bigger files, but compatible with current GMT scripts) or in binary format
# first receiver set (repeat these 6 lines and adjust nreceiversets accordingly)
-nrec = 1 # number of receivers
-xdeb = 1000. # first receiver x in meters
-zdeb = 1500. # first receiver z in meters
-xfin = 1000. # last receiver x in meters (ignored if onlyone receiver)
-zfin = 1500. # last receiver z in meters (ignored if onlyone receiver)
+nrec = 50 # number of receivers
+xdeb = 0.05 # first receiver x in meters
+zdeb = 0.2640 # first receiver z in meters
+xfin = 0.28 # last receiver x in meters (ignored if only one receiver)
+zfin = 0.2640 # last receiver z in meters (ignored if only one receiver)
enreg_surf_same_vertical = .false. # receivers inside the medium or at the surface
# display parameters
@@ -65,14 +66,14 @@ NSTEP_BETWEEN_OUTPUT_IMAGES = 100 # every how many time steps we
cutsnaps = 1. # minimum amplitude kept in % for the JPEG and PostScript snapshots; amplitudes below that are muted
#### for JPEG color images ####
output_color_image = .true. # output JPEG color image of the results every NSTEP_BETWEEN_OUTPUT_IMAGES time steps or not
-imagetype_JPEG = 10 # display 1=displ_Ux 2=displ_Uz 3=displ_norm 4=veloc_Vx 5=veloc_Vz 6=veloc_norm 7=accel_Ax 8=accel_Az 9=accel_norm 10=pressure
-factor_subsample_image = 1.0 # (double precision) factor to subsample color images output by the code (useful for very large models)
+imagetype_JPEG = 2 # display 1=displ_Ux 2=displ_Uz 3=displ_norm 4=veloc_Vx 5=veloc_Vz 6=veloc_norm 7=accel_Ax 8=accel_Az 9=accel_norm 10=pressure
+factor_subsample_image = 1.0 # (double precision) factor to subsample color images output by the code (useful for very large models)
POWER_DISPLAY_COLOR = 0.30d0 # non linear display to enhance small amplitudes in JPEG color images
DRAW_SOURCES_AND_RECEIVERS = .true. # display sources as orange crosses and receivers as green squares in JPEG images or not
DRAW_WATER_IN_BLUE = .true. # display acoustic layers as constant blue in JPEG images, because they likely correspond to water in the case of ocean acoustics or in the case of offshore oil industry experiments (if off, display them as greyscale, as for elastic or poroelastic elements, for instance for acoustic-only oil industry models of solid media)
USE_SNAPSHOT_NUMBER_IN_FILENAME = .false. # use snapshot number in the file name of JPEG color snapshots instead of the time step (for instance to create movies in an easier way later)
#### for PostScript snapshots ####
-output_postscript_snapshot = .false. # output Postscript snapshot of the results every NSTEP_BETWEEN_OUTPUT_IMAGES time steps or not
+output_postscript_snapshot = .true. # output Postscript snapshot of the results every NSTEP_BETWEEN_OUTPUT_IMAGES time steps or not
imagetype_postscript = 1 # display 1=displ vector 2=veloc vector 3=accel vector; small arrows are displayed for the vectors
meshvect = .true. # display mesh on PostScript plots or not
modelvect = .false. # display velocity model on PostScript plots or not
@@ -93,10 +94,10 @@ output_grid_ASCII = .false. # dump the grid in an ASCII tex
output_energy = .false. # compute and output total acoustic and elastic energy curves (slows down the code significantly)
# velocity and density models
-nbmodels = 4 # nb of different models
+nbmodels = 1 # nb of different models
# define models as
# I: (model_number 1 rho Vp Vs 0 0 QKappa Qmu 0 0 0 0 0 0) or
-# II: (model_number 2 rho c11 c13 c15 c33 c35 c55 0 0 0 0 0 0) or
+# II: (model_number 2 rho c11 c13 c15 c33 c35 c55 c12 c23 c25 0 0 0) or
# III: (model_number 3 rhos rhof phi c kxx kxz kzz Ks Kf Kfr etaf mufr Qmu).
#
# For istropic elastic/acoustic material use I and set Vs to zero to make a given model acoustic, for anisotropic elastic use II,
@@ -104,15 +105,21 @@ nbmodels = 4 # nb of different models
#
# For anisotropic media in the most general 2D case we have the following convention for the stress-strain relationship:
# ! implement anisotropy in 2D
-# sigma_xx = c11*dux_dx + c15*(duz_dx + dux_dz) + c13*duz_dz
-# sigma_zz = c13*dux_dx + c35*(duz_dx + dux_dz) + c33*duz_dz
-# sigma_xz = c15*dux_dx + c55*(duz_dx + dux_dz) + c35*duz_dz
+# sigma_xx = c11*dux_dx + c13*duz_dz + c15*(duz_dx + dux_dz)
+# sigma_zz = c13*dux_dx + c33*duz_dz + c35*(duz_dx + dux_dz)
+# sigma_xz = c15*dux_dx + c35*duz_dz + c55*(duz_dx + dux_dz)
+# sigma_yy is not equal to zero in the plane strain formulation
+# but is used only in post-processing if needed,
+# to compute pressure for display or seismogram recording purposes
+# sigma_yy = c12*dux_dx + c23*duz_dz + c25*(duz_dx + dux_dz)
# where the notations are for instance duz_dx = d(Uz) / dx.
#
-1 1 2700.d0 3000.d0 1732.051d0 0 0 9999 9999 0 0 0 0 0 0
-2 1 2500.d0 2700.d0 0 0 0 9999 9999 0 0 0 0 0 0
-3 1 2200.d0 2500.d0 1443.375d0 0 0 9999 9999 0 0 0 0 0 0
-4 1 2200.d0 2200.d0 1343.375d0 0 0 9999 9999 0 0 0 0 0 0
+# For anisotropic elastic media the last three parameters, c12 c23 c25, are used only when the user asks the code to compute pressure for display
+# or seismogram recording purposes. Thus, if you do not know these parameters for your anisotropic material and/or if you do not plan to display or record pressure you
+# can ignore them and set them to zero. When pressure is used these three parameters are needed because the code needs to compute sigma_yy,
+# which is not equal to zero in the plane strain formulation.
+#
+1 2 7100. 16.5d10 5.d10 0 6.2d10 0 3.96d10 0 0 0 0 0 0
# external mesh or not
read_external_mesh = .false.
@@ -151,12 +158,12 @@ tangential_detection_curve_file = ./DATA/courbe_eros_nodes # file containing the
# PARAMETERS FOR INTERNAL MESHING
# file containing interfaces for internal mesh
-interfacesfile = ./interfaces.dat
+interfacesfile = topoaniso.dat
# geometry of the model (origin lower-left corner = 0,0) and mesh description
-xmin = 0.d0 # abscissa of left side of the model
-xmax = 5000.d0 # abscissa of right side of the model
-nx = 100 # number of elements along X
+xmin = 0.0d0 # abscissa of left side of the model
+xmax = 0.33 # abscissa of right side of the model
+nx = 60 # number of elements along X
# absorbing boundary parameters (see absorbing_conditions above)
absorbbottom = .false.
@@ -165,7 +172,5 @@ absorbtop = .false.
absorbleft = .false.
# define the different regions of the model in the (nx,nz) spectral element mesh
-nbregions = 3 # nb of regions and model number for each
- 1 40 1 60 2
-41 60 1 60 1
-61 100 1 60 2
+nbregions = 1 # nb of regions and model number for each
+1 60 1 60 1
diff --git a/LuoYang_fluid_solid_kernel/SOURCE b/anisotropic_zinc_crystal/SOURCE
similarity index 87%
copy from LuoYang_fluid_solid_kernel/SOURCE
copy to anisotropic_zinc_crystal/SOURCE
index 528ca38..99ee461 100644
--- a/LuoYang_fluid_solid_kernel/SOURCE
+++ b/anisotropic_zinc_crystal/SOURCE
@@ -1,14 +1,14 @@
# The components of a moment tensor source must be given in N.m, not in dyne.cm as in the DATA/CMTSOLUTION source file of the 3D version of the code.
# source 1
source_surf = .false. # source inside the medium or at the surface
-xs = 4000. # source location x in meters
-zs = 1500. # source location z in meters
+xs = 0.165 # source location x in meters
+zs = 0.165 # source location z in meters
# source type: elastic force or acoustic pressure = 1 or moment tensor = 2;
# for a plane wave including converted and reflected waves at the free surface, P wave = 1, S wave = 2, Rayleigh wave = 3
# for a plane wave without converted nor reflected waves at the free surface, i.e. with the incident wave only, P wave = 4, S wave = 5
source_type = 1
time_function_type = 1 # Ricker = 1, first derivative = 2, Gaussian = 3, Dirac = 4, Heaviside = 5
-f0 = 10.0 # dominant source frequency (Hz) if not Dirac or Heaviside
+f0 = 170000. # dominant source frequency (Hz) if not Dirac or Heaviside
tshift = 0.0 # time shift when multi sources (if one source, must be zero)
anglesource = 0.0 # angle of the source (for a force only); for a plane wave, this is the incidence angle; for moment tensor sources this is unused
Mxx = 1. # Mxx component (for a moment tensor source only)
diff --git a/unused_older_examples_DATA_to_sort/interfaces_no_canyon.dat b/anisotropic_zinc_crystal/topoaniso.dat
similarity index 85%
copy from unused_older_examples_DATA_to_sort/interfaces_no_canyon.dat
copy to anisotropic_zinc_crystal/topoaniso.dat
index 0951555..048de0d 100644
--- a/unused_older_examples_DATA_to_sort/interfaces_no_canyon.dat
+++ b/anisotropic_zinc_crystal/topoaniso.dat
@@ -3,24 +3,24 @@
#
2
#
-# for each interface below, we give the number of points and then x,y for each point
+# for each interface below, we give the number of points and then x,z for each point
#
#
# interface number 1 (bottom of the mesh)
#
2
0 0
- 19 0
+ 0.33 0
#
# interface number 2 (topography, top of the mesh)
#
2
- 0 9
- 19 9
+ 0 0.33
+ 0.33 0.33
#
# for each layer, we give the number of spectral elements in the vertical direction
#
#
# layer number 1 (bottom layer)
#
- 45
+ 60
More information about the CIG-COMMITS
mailing list