Commit aef7f960 authored by Florian Wittkamp's avatar Florian Wittkamp

Squashed commit of the following:

Update: Merge test
parent 6761a267
This diff is collapsed.
/*-----------------------------------------------------------------------------------------
* Copyright (C) 2005 <name of author>
*
* This file is part of DENISE.
*
* DENISE is free software: you can redistribute it and/or modify
* it under the terms of the GNU General Public License as published by
* the Free Software Foundation, version 2.0 of the License only.
*
* DENISE is distributed in the hope that it will be useful,
* but WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
* GNU General Public License for more details.
*
* You should have received a copy of the GNU General Public License
* along with DENISE. See file COPYING and/or <http://www.gnu.org/licenses/gpl-2.0.html>.
-----------------------------------------------------------------------------------------*/
/* $Id: hh_el.c,v 2.3 2007/08/21 13:16:19 tbohlen Exp $*/
/*
* Model defined by flnode file. In this version a constant Q-model with Q_p=Q_s is used. It
* is defined via the relaxation frequencies and the TAU value given in the input file. The
* shear wave velocity model is defined by the values in file 'flnodes.Rheinstetten.start'
* whereas the P wave velocity model as well as the density are defined by the values in
* FLnode file 'flnodes.Rheinstetten'.
*/
#include "fd.h"
void model_elastic(float ** rho, float ** pi, float ** u){
/*--------------------------------------------------------------------------*/
/* extern variables */
extern int NX, NY, NXG, NYG, POS[3], MYID, L;
extern float DH, DT, *FL, TAU;
extern char MFILE[STRING_SIZE];
extern char INV_MODELFILE[STRING_SIZE];
/* local variables */
float muv, piv, vp, vs, rhov, ts, tp, *pts;
int i, j, ii, jj, l;
char modfile[STRING_SIZE];
FILE *flfile;
int nodes;
char cline[256];
float *fldepth, *flrho, *flvp, *flvs;
/**************************************************/
/* creation of shear wave velocity and TAU models */
/**************************************************/
/*read FL nodes from File*/
nodes=3;
fldepth=vector(1,nodes);
flvs=vector(1,nodes);
flrho=vector(1,nodes);
flvp=vector(1,nodes);
flfile=fopen("model_true/flnodes.toy_example.start","r");
if (flfile==NULL) err(" FL-file could not be opened !");
for (l=1;l<=nodes;l++){
fgets(cline,255,flfile);
if (cline[0]!='#'){
sscanf(cline,"%f%f%f%f",&fldepth[l], &flrho[l], &flvp[l], &flvs[l]);
}
else l=l-1;
}
if(MYID==0){
printf(" ------------------------------------------------------------------ \n\n");
printf(" Information of FL nodes: \n\n");
printf(" \t depth \t vs \n\n");
for (l=1;l<=nodes;l++){
printf(" \t %f \t %f\n\n",fldepth[l],flvs[l]);
}
printf(" ------------------------------------------------------------------ \n\n");
}
/* loop over global grid */
for (i=1;i<=NXG;i++){
for (l=1;l<nodes;l++){
if (fldepth[l]==fldepth[l+1]){
if ((i==1) && (MYID==0)){
printf("depth: %f m: double node\n",fldepth[l]);}}
else{
for (j=(int)(fldepth[l]/DH)+1;j<=(int)(fldepth[l+1]/DH);j++){
vs=0.0;
vs=(DH*(j-1)-fldepth[l])*(flvs[l+1]-flvs[l])/(fldepth[l+1]-fldepth[l])+flvs[l];
vs=vs*1000.0;
muv=vs;
ts=TAU;
tp=TAU;
/* only the PE which belongs to the current global gridpoint
is saving model parameters in his local arrays */
if ((POS[1]==((i-1)/NX)) &&
(POS[2]==((j-1)/NY))){
ii=i-POS[1]*NX;
jj=j-POS[2]*NY;
u[jj][ii]=muv;
}
}
}
}
for (j=(int)(fldepth[nodes]/DH)+1;j<=NYG;j++){
vs=0.0;
vs=flvs[nodes]*1000.0;
muv=vs;
/* only the PE which belongs to the current global gridpoint
is saving model parameters in his local arrays */
if ((POS[1]==((i-1)/NX)) &&
(POS[2]==((j-1)/NY))){
ii=i-POS[1]*NX;
jj=j-POS[2]*NY;
u[jj][ii]=muv;
}
}
}
free_vector(fldepth,1,nodes);
free_vector(flrho,1,nodes);
free_vector(flvp,1,nodes);
free_vector(flvs,1,nodes);
/**************************************************/
/* creation of P wave velocity and density models */
/**************************************************/
/*read FL nodes from File*/
nodes=7;
fldepth=vector(1,nodes);
flvs=vector(1,nodes);
flrho=vector(1,nodes);
flvp=vector(1,nodes);
flfile=fopen("model_true/flnodes.toy_example","r");
if (flfile==NULL) err(" FL-file could not be opened !");
for (l=1;l<=nodes;l++){
fgets(cline,255,flfile);
if (cline[0]!='#'){
sscanf(cline,"%f%f%f%f",&fldepth[l], &flrho[l], &flvp[l], &flvs[l]);
}
else l=l-1;
}
if(MYID==0){
printf(" ------------------------------------------------------------------ \n\n");
printf(" Information of FL nodes: \n\n");
printf(" \t depth \t vp \t rho \n\n");
for (l=1;l<=nodes;l++){
printf(" \t %f \t %f \t %f \n\n",fldepth[l],flvp[l],flrho[l]);
}
printf(" ------------------------------------------------------------------ \n\n");
}
/* loop over global grid */
for (i=1;i<=NXG;i++){
for (l=1;l<nodes;l++){
if (fldepth[l]==fldepth[l+1]){
if ((i==1) && (MYID==0)){
printf("depth: %f m: double node\n",fldepth[l]);}}
else{
for (j=(int)(fldepth[l]/DH)+1;j<=(int)(fldepth[l+1]/DH);j++){
vp=0.0; rhov=0.0;
vp=(DH*(j-1)-fldepth[l])*(flvp[l+1]-flvp[l])/(fldepth[l+1]-fldepth[l])+flvp[l];
vp=vp*1000.0;
rhov=(DH*(j-1)-fldepth[l])*(flrho[l+1]-flrho[l])/(fldepth[l+1]-fldepth[l])+flrho[l];
rhov=rhov*1000.0;
piv=vp;
/* only the PE which belongs to the current global gridpoint
is saving model parameters in his local arrays */
if ((POS[1]==((i-1)/NX)) &&
(POS[2]==((j-1)/NY))){
ii=i-POS[1]*NX;
jj=j-POS[2]*NY;
rho[jj][ii]=rhov;
pi[jj][ii]=piv;
}
}
}
}
for (j=(int)(fldepth[nodes]/DH)+1;j<=NYG;j++){
vp=0.0; rhov=0.0;
vp=flvp[nodes]*1000.0; rhov=flrho[nodes]*1000.0;
piv=vp;
/* only the PE which belongs to the current global gridpoint
is saving model parameters in his local arrays */
if ((POS[1]==((i-1)/NX)) &&
(POS[2]==((j-1)/NY))){
ii=i-POS[1]*NX;
jj=j-POS[2]*NY;
rho[jj][ii]=rhov;
pi[jj][ii]=piv;
}
}
}
free_vector(fldepth,1,nodes);
free_vector(flrho,1,nodes);
free_vector(flvp,1,nodes);
free_vector(flvs,1,nodes);
sprintf(modfile,"%s_rho_it0.bin",INV_MODELFILE);
writemod(modfile,rho,3);
MPI_Barrier(MPI_COMM_WORLD);
if (MYID==0) mergemod(modfile,3);
MPI_Barrier(MPI_COMM_WORLD);
sprintf(modfile,"%s_rho_it0.bin.%i%i",INV_MODELFILE,POS[1],POS[2]);
remove(modfile);
sprintf(modfile,"%s_vs_it0.bin",INV_MODELFILE);
writemod(modfile,u,3);
MPI_Barrier(MPI_COMM_WORLD);
if (MYID==0) mergemod(modfile,3);
MPI_Barrier(MPI_COMM_WORLD);
sprintf(modfile,"%s_vs_it0.bin.%i%i",INV_MODELFILE,POS[1],POS[2]);
remove(modfile);
sprintf(modfile,"%s_vp_it0.bin",INV_MODELFILE);
writemod(modfile,pi,3);
MPI_Barrier(MPI_COMM_WORLD);
if (MYID==0) mergemod(modfile,3);
MPI_Barrier(MPI_COMM_WORLD);
sprintf(modfile,"%s_vp_it0.bin.%i%i",INV_MODELFILE,POS[1],POS[2]);
remove(modfile);
free_vector(pts,1,L);
}
/*-----------------------------------------------------------------------------------------
* Copyright (C) 2005 <name of author>
*
* This file is part of DENISE.
*
* DENISE is free software: you can redistribute it and/or modify
* it under the terms of the GNU General Public License as published by
* the Free Software Foundation, version 2.0 of the License only.
*
* DENISE is distributed in the hope that it will be useful,
* but WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
* GNU General Public License for more details.
*
* You should have received a copy of the GNU General Public License
* along with DENISE. See file COPYING and/or <http://www.gnu.org/licenses/gpl-2.0.html>.
-----------------------------------------------------------------------------------------*/
/*
* Model defined by flnode file.
*/
#include "fd.h"
void model_elastic(float ** rho, float ** pi, float ** u){
/*--------------------------------------------------------------------------*/
/* extern variables */
extern int NX, NY, NXG, NYG, POS[3], MYID, L;
extern float DH, DT, *FL, TAU;
extern char MFILE[STRING_SIZE];
extern char INV_MODELFILE[STRING_SIZE];
/* local variables */
float muv, piv, vp, vs, rhov, ts, tp, *pts;
int i, j, ii, jj, l;
char modfile[STRING_SIZE];
FILE *flfile;
int nodes;
char cline[256];
float *fldepth, *flrho, *flvp, *flvs;
nodes=7;
fldepth=vector(1,nodes);
flrho=vector(1,nodes);
flvp=vector(1,nodes);
flvs=vector(1,nodes);
/*read FL nodes from File*/
flfile=fopen("model_true/flnodes.toy_example","r");
if (flfile==NULL) err(" FL-file could not be opened !");
for (l=1;l<=nodes;l++){
fgets(cline,255,flfile);
if (cline[0]!='#'){
sscanf(cline,"%f%f%f%f",&fldepth[l], &flrho[l], &flvp[l], &flvs[l]);
}
else l=l-1;
}
if(MYID==0){
printf(" ------------------------------------------------------------------ \n\n");
printf(" Information of FL nodes: \n\n");
printf(" \t depth \t rho \t vp \t vs \n\n");
for (l=1;l<=nodes;l++){
printf(" \t %f \t %f \t %f \t %f\n\n",fldepth[l],flrho[l],flvp[l],flvs[l]);
}
printf(" ------------------------------------------------------------------ \n\n");
}
/*-----------------------------------------------------------------------*/
/*-----------------------------------------------------------------------*/
/* loop over global grid */
for (i=1;i<=NXG;i++){
for (l=1;l<nodes;l++){
if (fldepth[l]==fldepth[l+1]){
if ((i==1) && (MYID==0)){
printf("depth: %f m: double node\n",fldepth[l]);}}
else{
for (j=(int)(fldepth[l]/DH)+1;j<=(int)(fldepth[l+1]/DH);j++){
vp=0.0; vs=0.0; rhov=0.0;
vp=(DH*(j-1)-fldepth[l])*(flvp[l+1]-flvp[l])/(fldepth[l+1]-fldepth[l])+flvp[l];
vp=vp*1000.0;
vs=(DH*(j-1)-fldepth[l])*(flvs[l+1]-flvs[l])/(fldepth[l+1]-fldepth[l])+flvs[l];
vs=vs*1000.0;
rhov=(DH*(j-1)-fldepth[l])*(flrho[l+1]-flrho[l])/(fldepth[l+1]-fldepth[l])+flrho[l];
rhov=rhov*1000.0;
muv=vs;
piv=vp;
ts=TAU;
tp=TAU;
/* only the PE which belongs to the current global gridpoint
is saving model parameters in his local arrays */
if ((POS[1]==((i-1)/NX)) &&
(POS[2]==((j-1)/NY))){
ii=i-POS[1]*NX;
jj=j-POS[2]*NY;
u[jj][ii]=muv;
rho[jj][ii]=rhov;
pi[jj][ii]=piv;
}
}
}
}
for (j=(int)(fldepth[nodes]/DH)+1;j<=NYG;j++){
vp=0.0; vs=0.0; rhov=0.0;
vp=flvp[nodes]*1000.0; vs=flvs[nodes]*1000.0; rhov=flrho[nodes]*1000.0;
muv=vs;
piv=vp;
ts=TAU;
tp=TAU;
/* only the PE which belongs to the current global gridpoint
is saving model parameters in his local arrays */
if ((POS[1]==((i-1)/NX)) &&
(POS[2]==((j-1)/NY))){
ii=i-POS[1]*NX;
jj=j-POS[2]*NY;
u[jj][ii]=muv;
rho[jj][ii]=rhov;
pi[jj][ii]=piv;
}
}
}
sprintf(modfile,"%s_rho_it0.bin",INV_MODELFILE);
writemod(modfile,rho,3);
MPI_Barrier(MPI_COMM_WORLD);
if (MYID==0) mergemod(modfile,3);
MPI_Barrier(MPI_COMM_WORLD);
sprintf(modfile,"%s_rho_it0.bin.%i%i",INV_MODELFILE,POS[1],POS[2]);
remove(modfile);
sprintf(modfile,"%s_vs_it0.bin",INV_MODELFILE);
writemod(modfile,u,3);
MPI_Barrier(MPI_COMM_WORLD);
if (MYID==0) mergemod(modfile,3);
MPI_Barrier(MPI_COMM_WORLD);
sprintf(modfile,"%s_vs_it0.bin.%i%i",INV_MODELFILE,POS[1],POS[2]);
remove(modfile);
sprintf(modfile,"%s_vp_it0.bin",INV_MODELFILE);
writemod(modfile,pi,3);
MPI_Barrier(MPI_COMM_WORLD);
if (MYID==0) mergemod(modfile,3);
MPI_Barrier(MPI_COMM_WORLD);
sprintf(modfile,"%s_vp_it0.bin.%i%i",INV_MODELFILE,POS[1],POS[2]);
remove(modfile);
free_vector(fldepth,1,nodes);
free_vector(flrho,1,nodes);
free_vector(flvp,1,nodes);
free_vector(flvs,1,nodes);
free_vector(pts,1,L);
}
......@@ -72,4 +72,7 @@
"General inversion parameters" : "comment",
"INVMAT1" : "1",
"INVMAT" : "10",
"Verbose mode" : "comment",
"VERBOSE" : "0",
}
......@@ -6,6 +6,7 @@
#
{
"Domain Decomposition" : "comment",
"NPROCX" : "4",
"NPROCY" : "2",
......@@ -96,4 +97,7 @@
"General inversion parameters" : "comment",
"INVMAT1" : "1",
"INVMAT" : "10",
"Verbose mode" : "comment",
"VERBOSE" : "0",
}
......@@ -104,4 +104,6 @@
"Termination of the programmme" : "comment",
"PRO" : "0.01",
"Verbose mode" : "comment",
"VERBOSE" : "0",
}
......@@ -225,4 +225,7 @@
"TWLENGTH_PLUS" : "0.01",
"TWLENGTH_MINUS" : "0.01",
"GAMMA" : "100000",
"Verbose mode" : "comment",
"VERBOSE" : "0",
}
......@@ -80,4 +80,6 @@
"INVMAT1" : "1",
"INVMAT" : "10",
"Verbose mode" : "comment",
"VERBOSE" : "0",
}
#-----------------------------------------------------------------
# JSON PARAMETER FILE FOR DENISE
#-----------------------------------------------------------------
# description:
# description/name of the model: toy_example_true.c with flnodes.toy_example
#
{
"Domain Decomposition" : "comment",
"NPROCX" : "4",
"NPROCY" : "1",
"FD order" : "comment",
"FDORDER" : "4",
"MAXRELERROR" : "0",
"2-D Grid" : "comment",
"NX" : "252",
"NY" : "75",
"DH" : "0.2",
"Time Stepping" : "comment",
"TIME" : "0.6",
"DT" : "5.0e-05",
"Source" : "comment",
"QUELLART" : "4",
"QUELLART_SH" : "4",
"QUELLART values: ricker=1;fumue=2;from_SIGNAL_FILE=3;SIN**3=4;Gaussian_deriv=5;Spike=6" : "comment",
"SIGNAL_FILE" : "./ormsby.dat",
"SIGNAL_FILE_SH" : "./ormsby.dat",
"QUELLTYP" : "3",
"QUELLTYP values (point_source) P-SV: explosive=1;force_in_x=2;force_in_y=3;force_in_z=5;rotated_force=4" : "comment",
"SRCREC" : "1",
"SRCREC values : read source positions from SOURCE_FILE=1, PLANE_WAVE=2" : "comment",
"SOURCE_FILE" : "./source/src_toy_example.dat",
"RUN_MULTIPLE_SHOTS" : "1",
"Acoustic Computation" : "comment",
"ACOUSTIC" : "0",
"Wavetype Computation" : "comment",
"WAVETYPE" : "3",
"Verbose" : "comment",
"VERBOSE" : "0",
"Model" : "comment",
"READMOD" : "0",
"MFILE" : "model/mod_toy_example_true",
"Free Surface" : "comment",
"FREE_SURF" : "1",
"PML Boundary" : "comment",
"FW" : "10",
"DAMPING" : "1700.0",
"FPML" : "31.25",
"BOUNDARY" : "0",
"npower" : "4.0",
"k_max_PML" : "1.0",
"Receiver" : "comment",
"SEISMO" : "5",
"READREC" : "1",
"REC_FILE" : "./receiver/rec_toy_example.dat",
"REFRECX, REFRECY" : "0.0 , 0.0",
"Seismograms" : "comment",
"NDT" : "1",
"SEIS_FORMAT" : "1",
"SEIS_FILE_P" : "su/measured_data/toy_example_p.su",
"SEIS_FILE_VX" : "su/measured_data/toy_example_vx.su",
"SEIS_FILE_VY" : "su/measured_data/toy_example_vy.su",
"SEIS_FILE_VZ" : "su/measured_data/toy_example_vz.su",
"Q-approximation" : "comment",
"L" : "0",
"FL1" : "0.5211",
"FL2" : "7.6660",
"FL3" : "72.6774",
"TAU" : "0.0966",
"General inversion parameters" : "comment",
"ITERMAX" : "1",
"INVMAT1" : "1",
"INVMAT" : "10",
}
......@@ -82,7 +82,12 @@
"INVMAT" : "0",
"QUELLTYPB" : "1",
"MISFIT_LOG_FILE" : "LOG_toy_example.dat",
"Workflow" : "comment",
"USE_WORKFLOW" : "1",
"FILE_WORKFLOW" : "workflow.txt",
"RESTART_WORKFLOW" : "0",
"Inversion for density" : "comment",
"INV_RHO_ITER" : "400",
"INV_VP_ITER" : "400",
......@@ -145,4 +150,7 @@
"FC_END" : "70.0",
"FC_INCR" : "10.0",
"ORDER" : "4",
"Verbose mode" : "comment",
"VERBOSE" : "0",
}
#-----------------------------------------------------------------
# JSON PARAMETER FILE FOR DENISE
#-----------------------------------------------------------------
# description:
# description/name of the model: genmod/toy_example_start.c with flnodes.toy_example.start
#
{
"Domain Decomposition" : "comment",
"NPROCX" : "4",
"NPROCY" : "1",
"FD order" : "comment",
"FDORDER" : "4",
"MAXRELERROR" : "0",
"2-D Grid" : "comment",
"NX" : "252",
"NY" : "75",
"DH" : "0.2",
"Time Stepping" : "comment",
"TIME" : "0.6",
"DT" : "5.0e-05",
"Source" : "comment",
"QUELLART" : "4",
"QUELLART_SH" : "4",
"QUELLART values: ricker=1;fumue=2;from_SIGNAL_FILE=3;SIN**3=4;Gaussian_deriv=5;Spike=6" : "comment",
"SIGNAL_FILE" : "./ormsby.dat",
"SIGNAL_FILE_SH" : "./ormsby.dat",
"QUELLTYP" : "3",
"QUELLTYP values (point_source) P-SV: explosive=1;force_in_x=2;force_in_y=3;rotated_force=4" : "comment",
"SRCREC" : "1",
"SRCREC values : read source positions from SOURCE_FILE=1, PLANE_WAVE=2" : "comment",
"SOURCE_FILE" : "./source/src_toy_example.dat",
"RUN_MULTIPLE_SHOTS" : "1",
"Acoustic Computation" : "comment",
"ACOUSTIC" : "0",
"Wavetype Computation" : "comment",
"WAVETYPE" : "2",
"Verbose" : "comment",
"VERBOSE" : "0",