[CIG-SEISMO] gmsh and specfem2D
giuseppe.digiulio at ingv.it
giuseppe.digiulio at ingv.it
Mon Sep 24 09:23:00 PDT 2012
I thank Paul Cristini to have fixed my error on *geo file created by
means of gmsh.
Now I have two doubts i) on the element sizes of my mesh and ii) on
the use of specfem2D with a planar incident vertically waves.
I built my external mesh by the command
gmsh 3UFscale_allquad.geo -2 -clmin 100 -clmax 200 -order 2
My *geo file is composed of three uniform parallel layers (with two
bottom layers having the same physical characteristics).
The total dimension of my working space is a 1000x1000 box.
I imagine with the above gmsh command to fix the minimun and maximum
length of my quad elements to 100 and 200, respectively. But my mesh
shows 20 elemets for each boundary side (see file
3UFscale_allquad.msh.jpg); it looks the dimension of quad elements is
shorter than 100.
The second doubts is on the use of a planar incident vertically waves
with specfem2D in conjunction with my external mesh (in other words
the paragraph how to set plane waves as initial conditions in the
specfem2d manual).
My Par_file assumes PSV waves, and the source_type is a SV waves. I am
looking to fix the plane wavefield at the interface between second
and third layer (at a depth of 500 m following my mesh). The planar
wavefield should be vertically incident. Looking at the initial
vector*ps files, it seems that I have two waves at the bottom and at
the top of my models, and not inside at the depth of 500 m as I would
like.
I attach my trial
Thanks in advance
Giuseppe Di Giulio
----------------------------------------------------------------
Il contenuto di questa e-mail e' rivolto unicamente alle persone cui
e' indirizzato, e puo'
contenere informazioni la cui riservatezza e' tutelata.
E' proibita la copia, la divulgazione o l'uso di questo messaggio o
dell'informazione ivi
contenuta da chiunque altro che non sia il destinatario indicato.
Se avete ricevuto questa e-mail per errore, vogliate cortesemente
comunicarlo immediatamente per telefono, fax o e-mail. Grazie.
This e-mail is intended only for person or entity to which is
addressed and may contain
information that is privileged, confidential or otherwise protected
from disclosure.
Copying, dissemination or use of this e-mail or the information herein
by anyone other than the intended recipient is prohibited.
If you have received this
e-mail by mistake, please notify us immediately by telephone, fax or e-mail.
Thank you.
-------------- next part --------------
lc = 1e-1;
Point(1) = {0, 0, 0, lc};
Point(2) = {1000, 0, 0, lc};
Point(3) = {0, 1000, 0, lc};
Point(4) = {1000, 1000, 0, lc};
Point(5) = {0, 500, 0, lc};
Point(6) = {1000, 500, 0, lc};
Point(7) = {0, 800, 0, lc};
Point(8) = {1000, 800, 0, lc};
Line(1) = {1, 2};
Line(2) = {2, 6};
Line(3) = {6, 5};
Line(4) = {5, 1};
Line(5) = {5, 7};
Line(6) = {7, 8};
Line(7) = {8, 6};
Line(8) = {8, 4};
Line(9) = {4, 3};
Line(10) = {3, 7};
Line Loop(11) = {1, 2, 3, 4};
Plane Surface(12) = {11};
Line Loop(13) = {-3, -5, -6, -7};
Plane Surface(14) = {13};
Line Loop(15) = {6, 8, 9, 10};
Plane Surface(16) = {15};
Recombine Surface{12,14,16};
Mesh.SubdivisionAlgorithm=1;
Mesh.RandomFactor=1e-6;
Physical Line("Top") = {9};
Physical Line("Left") = {10,-5,4};
Physical Line("Bottom") = {1};
Physical Line("Right") = {2,-7,8};
Physical Surface("M1") = {16};
Physical Surface("M2") = {14};
Physical Surface("M3") = {12};
-------------- next part --------------
A non-text attachment was scrubbed...
Name: 3UFscale_allquad.msh
Type: model/mesh
Size: 99396 bytes
Desc: not available
Url : http://geodynamics.org/pipermail/cig-seismo/attachments/20120924/1123de70/attachment-0001.msh
-------------- next part --------------
# title of job, and file that contains interface data r19201
title = 3UF Three Uniform layers
# 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 = .false. # save the last frame, needed for adjoint simulation
# parameters concerning partitioning
nproc = 1 # number of processes
partitioning_method = 3 # ascending order = 1, Scotch = 3
PERFORM_CUTHILL_MCKEE = .true. # perform inverse Cuthill-McKee (1969) optimization/permutation for mesh numbering
ngnod = 4 # number of control nodes per element (4 or 9)
initialfield = .true. # 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
READ_EXTERNAL_SEP_FILE = .false. # Read external SEP file from DATA/model_velocity.dat_input, or use routine
ATTENUATION_VISCOELASTIC_SOLID = .false. # turn attenuation (viscoelasticity) on or off for non-poroelastic solid parts of the model
ATTENUATION_PORO_FLUID_PART = .false. # turn viscous attenuation on or off for the fluid part of poroelastic parts of the model
Q0 = 1 # quality factor for viscous attenuation
freq0 = 10 # frequency for viscous attenuation
p_sv = .true. # set the type of calculation (P-SV or SH/membrane waves)
# time step parameters
nt = 80000 # total number of time steps
deltat = 1.d-4 # duration of a time step
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 4th-order 4-stage Runge-Kutta
# source parameters
NSOURCES = 1 # number of sources [source info read in CMTSOLUTION file]
force_normal_to_surface = .false. # angleforce normal to surface (external mesh and curve file needed)
# constants for attenuation
N_SLS = 2 # number of standard linear solids for attenuation
f0_attenuation = 5.196152422706633 # (Hz) relevant only if source is a Dirac or a Heaviside, else it is f0
# receiver set parameters for seismograms
seismotype = 1 # record 1=displ 2=veloc 3=accel 4=pressure
generate_STATIONS = .true. # creates a STATION file in ./DATA
nreceiversets = 1 # number of receiver sets
anglerec = 0.d0 # angle to rotate components at receivers
rec_normal_to_surface = .false. # replace anglerec with the normal to the surface (external mesh and curve file needed, anglerec above ignored)
SU_FORMAT = .false. # output seismograms in Seismic Unix format (adjoint traces will be read in the same format)
# first receiver set
nrec = 10 # number of receivers
xdeb = 10. # first receiver x in meters
zdeb = 1000 # first receiver z in meters
xfin = 900 # last receiver x in meters (ignored if onlyone receiver)
zfin = 1000 # last receiver z in meters (ignored if onlyone receiver)
enreg_surf_same_vertical = .false. # receivers inside the medium or at the surface
# display parameters
NTSTEP_BETWEEN_OUTPUT_INFO = 200 # display frequency in time steps
output_postscript_snapshot = .true. # output Postscript snapshot of the results
output_color_image = .true. # output color image of the results
imagetype = 1 # display 1=displ 2=veloc 3=accel 4=pressure
cutsnaps = 1. # minimum amplitude in % for snapshots
meshvect = .true. # display mesh on vector plots or not
modelvect = .true. # display velocity model on vector plots
boundvect = .true. # display boundary conditions on plots
interpol = .false. # interpolation of the display or not
pointsdisp = 6 # points for interpolation of display (set to 1 for lower-left corner only)
subsamp_postscript = 8 # subsampling of color snapshots
factor_subsample_image = 3 # 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 color images
DRAW_WATER_CONSTANT_BLUE_IN_JPG = .true. # display acoustic layers as constant blue in JPEG images, because they likely correspond to water
sizemax_arrows = 1.d0 # maximum size of arrows on vector plots in cm
US_LETTER = .false. # US letter paper or European A4
USE_SNAPSHOT_NUMBER_IN_FILENAME = .false. # use snapshot number in the file name of JPEG color snapshots instead of the time step
gnuplot = .false. # generate a GNUPLOT file for the grid
output_grid = .false. # save the grid in a text file or not
output_energy = .false. # compute and output acoustic and elastic energy (slows down the code significantly)
output_wavefield_snapshot = .false. # output Ux,Uy,Uz text file for each output time (big files)
# velocity and density models
nbmodels = 3 # nb of different models
# define models as I: (model_number,1,rho,Vp,Vs,0,0,QKappa Qmu) or II: (model_number,2,rho,c11,c13,c15,c33,c35,c55) 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,
# and for isotropic poroelastic material use III. The mesh can contain acoustic, elastic, & poroelastic models simultaneously
1 1 2000.d0 1000.d0 500.0 0 0 100. 50. 0 0 0 0 0 0
2 1 2000.d0 1000.d0 500.0 0 0 100. 50. 0 0 0 0 0 0
3 1 2000.d0 3000.d0 1500.0 0 0 300. 150. 0 0 0 0 0 0
#3 3 2650.d0 1000.0d0 0.4 1.25 1d-10 0.0 1d-10 3.6d10 2.25d9 2d9 0.0d-4 3.204d9 10.d0
#2 1 2500.d0 5000.d0 2500.0d0 0 0 10.d0 10.d0 0 0 0 0 0 0
#3 1 2200.d0 2500.d0 1443.375d0 0 0 10.d0 10.d0 0 0 0 0 0 0
#4 3 2200.d0 786.3d0 0.4 2.0 1d-11 0.0 1d-11 5.341d9 2d9 3d9 0.0d-4 3.204d9 10.d0
#5 2 2500.d0 169.d9 122.d9 0.d0 169.d9 0.d0 75.3d9 0 0 0 0 0 0
#4 1 2200.d0 2200.d0 1343.375d0 0 0 10.d0 10.d0 0 0 0 0 0 0
# external mesh or not
read_external_mesh = .true.
# absorbing boundary active or not
absorbing_conditions = .true.
# for horizontal periodic conditions: detect common points between left and right edges
ADD_PERIODIC_CONDITIONS = .false.
# horizontal periodicity distance for periodic conditions
PERIODIC_horiz_dist = 0.3597d0
# grid point detection tolerance for periodic conditions
PERIODIC_DETECT_TOL = 3.3334d-6
#-----------------------------------------------------------------------------
# PARAMETERS FOR EXTERNAL MESHING
# data concerning mesh, when generated using third-party app (more info in README)
# (see also absorbing_conditions above)
mesh_file = Mesh_3UFscale_allquad # file containing the mesh
nodes_coords_file = Nodes_3UFscale_allquad # file containing the nodes coordinates
materials_file = Material_3UFscale_allquad # file containing the material number for each element
free_surface_file = Surf_free_3UFscale_allquad # file containing the free surface
absorbing_surface_file = Surf_abs_3UFscale_allquad # file containing the absorbing surface
tangential_detection_curve_file = ./DATA/courbe_eros_nodes # file containing the curve delimiting the velocity model
#-----------------------------------------------------------------------------
# PARAMETERS FOR INTERNAL MESHING
# file containing interfaces for internal mesh
interfacesfile = Interface_flat_ASM_DGA_119_62kHz.dat
# geometry of the model (origin lower-left corner = 0,0) and mesh description
xmin = -1.d0 # abscissa of left side of the model
xmax = 4.d0 # abscissa of right side of the model# file containing interfaces for internal mesh
nx = 835 # number of elements along X
# absorbing boundaries parameters
absorbbottom = .true.
absorbright = .true.
absorbtop = .true.
absorbleft = .true.
# define the different regions of the model in the (nx,nz) spectral element mesh
nbregions = 2 # nb of regions and model number for each
1 835 1 71 1
1 835 72 96 2
#1 80 21 40 4
#1 80 41 60 3
#60 80 21 40 4
#30 40 50 60 2
#35 40 50 60 5
-------------- next part --------------
#source 1. 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_surf = .true. # source inside the medium or at the surface
xs = 0. # source location x in meters
zs = 500 # source location z in meters
source_type = 2 # elastic force or acoustic pressure = 1 or moment tensor = 2
time_function_type = 4 # Ricker = 1, first derivative = 2, Gaussian = 3, Dirac = 4, Heaviside = 5
f0 = 0.25d0 # dominant source frequency (Hz) if not Dirac or Heaviside
t0 = 0.0 # time shift when multi sources (if one source, must be zero)
angleforce = 0.0 # angle of the source (for a force only)
Mxx = 1. # Mxx component (for a moment tensor source only)
Mzz = 1. # Mzz component (for a moment tensor source only)
Mxz = 0. # Mxz component (for a moment tensor source only)
factor = 1.d10 # amplification factor
-------------- next part --------------
A non-text attachment was scrubbed...
Name: my_example.tar.gz
Type: application/x-gzip
Size: 1723993 bytes
Desc: not available
Url : http://geodynamics.org/pipermail/cig-seismo/attachments/20120924/1123de70/attachment-0001.bin
More information about the CIG-SEISMO
mailing list