[cig-commits] r8433 - in seismo/2D/SPECFEM2D/trunk: MAILLE90 MAILLE90/old_topo_files SPECFEM90

walter at geodynamics.org walter at geodynamics.org
Fri Dec 7 15:45:56 PST 2007


Author: walter
Date: 2007-12-07 15:45:56 -0800 (Fri, 07 Dec 2007)
New Revision: 8433

Added:
   seismo/2D/SPECFEM2D/trunk/MAILLE90/old_topo_files/
   seismo/2D/SPECFEM2D/trunk/MAILLE90/old_topo_files/interf_paco.dat
   seismo/2D/SPECFEM2D/trunk/MAILLE90/old_topo_files/topo_cours_M2_UPPA_curved.dat
   seismo/2D/SPECFEM2D/trunk/MAILLE90/old_topo_files/topo_rouen_aniso.dat
   seismo/2D/SPECFEM2D/trunk/MAILLE90/old_topo_files/topoarticle.dat
Removed:
   seismo/2D/SPECFEM2D/trunk/MAILLE90/interf_paco.dat
   seismo/2D/SPECFEM2D/trunk/MAILLE90/profilx.dat
   seismo/2D/SPECFEM2D/trunk/MAILLE90/profily.dat
   seismo/2D/SPECFEM2D/trunk/MAILLE90/topo_cours_M2_UPPA_curved.dat
   seismo/2D/SPECFEM2D/trunk/MAILLE90/topo_rouen_aniso.dat
   seismo/2D/SPECFEM2D/trunk/MAILLE90/topoarticle.dat
Modified:
   seismo/2D/SPECFEM2D/trunk/MAILLE90/Makefile
   seismo/2D/SPECFEM2D/trunk/MAILLE90/Par
   seismo/2D/SPECFEM2D/trunk/MAILLE90/Par_cours_M2_UPPA
   seismo/2D/SPECFEM2D/trunk/MAILLE90/meshfem2D.f90
   seismo/2D/SPECFEM2D/trunk/SPECFEM90/Makefile
   seismo/2D/SPECFEM2D/trunk/SPECFEM90/checkgrid.f90
   seismo/2D/SPECFEM2D/trunk/SPECFEM90/compute_gradient_attenuation.f90
   seismo/2D/SPECFEM2D/trunk/SPECFEM90/createnum_fast.f90
   seismo/2D/SPECFEM2D/trunk/SPECFEM90/createnum_slow.f90
   seismo/2D/SPECFEM2D/trunk/SPECFEM90/datim.f90
   seismo/2D/SPECFEM2D/trunk/SPECFEM90/defarrays.f90
   seismo/2D/SPECFEM2D/trunk/SPECFEM90/define_derivative_matrices.f90
   seismo/2D/SPECFEM2D/trunk/SPECFEM90/gll_library.f90
   seismo/2D/SPECFEM2D/trunk/SPECFEM90/gmat01.f90
   seismo/2D/SPECFEM2D/trunk/SPECFEM90/lagrange_poly.f90
   seismo/2D/SPECFEM2D/trunk/SPECFEM90/plotgll.f90
   seismo/2D/SPECFEM2D/trunk/SPECFEM90/plotpost.f90
   seismo/2D/SPECFEM2D/trunk/SPECFEM90/positrec.f90
   seismo/2D/SPECFEM2D/trunk/SPECFEM90/positsource.f90
   seismo/2D/SPECFEM2D/trunk/SPECFEM90/q49shape.f90
   seismo/2D/SPECFEM2D/trunk/SPECFEM90/q49spec.f90
   seismo/2D/SPECFEM2D/trunk/SPECFEM90/specfem2D.f90
   seismo/2D/SPECFEM2D/trunk/SPECFEM90/write_seismograms.f90
Log:
made 2D mesher more general, can now handle any number of curved layers.
also cleaned a few details in the solver.


Modified: seismo/2D/SPECFEM2D/trunk/MAILLE90/Makefile
===================================================================
--- seismo/2D/SPECFEM2D/trunk/MAILLE90/Makefile	2004-12-14 15:22:01 UTC (rev 8432)
+++ seismo/2D/SPECFEM2D/trunk/MAILLE90/Makefile	2007-12-07 23:45:56 UTC (rev 8433)
@@ -1,7 +1,7 @@
 #
-# Makefile for MESHFEM2D version 5.0
+# Makefile for MESHFEM2D version 5.1
 #
-# Dimitri Komatitsch, Universite de Pau et des Pays de l'Adour, May 2004
+# Dimitri Komatitsch, Universite de Pau et des Pays de l'Adour, December 2004
 # 
 SHELL=/bin/sh
 
@@ -13,8 +13,8 @@
 
 # Intel Linux
 F90 = ifort
-#FLAGS=-O0 -implicitnone -warn stderrors -warn truncated_source -warn argument_checking -warn unused -warn declarations -std95 -check bounds -C
-FLAGS=-O2 -implicitnone -warn stderrors -warn truncated_source -warn argument_checking -warn unused -warn declarations -std95 -check nobounds
+FLAGS=-O0 -implicitnone -warn stderrors -warn truncated_source -warn argument_checking -warn unused -warn declarations -std95 -check bounds -CB
+#FLAGS=-O2 -implicitnone -warn stderrors -warn truncated_source -warn argument_checking -warn unused -warn declarations -std95 -check nobounds
 
 #
 # g95 (free f95 compiler from http://www.g95.org, still under development, but works)

Modified: seismo/2D/SPECFEM2D/trunk/MAILLE90/Par
===================================================================
--- seismo/2D/SPECFEM2D/trunk/MAILLE90/Par	2004-12-14 15:22:01 UTC (rev 8432)
+++ seismo/2D/SPECFEM2D/trunk/MAILLE90/Par	2007-12-07 23:45:56 UTC (rev 8433)
@@ -6,20 +6,17 @@
 # ----------------------------------------------------------------
 <-                              ->
 #
-# File names and path for different outputs
+# title of job, and file that contains interface data
 #
 title                           = Test pour cours M2 UPPA
-topofile                        = topo_cours_M2_UPPA_curved.dat
-interffile                      = none
+interfacesfile                  = interfaces_cours_M2_UPPA_curved.dat
 #
 # geometry of the model (origin lower-left corner = 0,0) and mesh description
 #
-xmin                            = 0.0d0          ! abscissa of left side of the model
+xmin                            = 0.d0           ! abscissa of left side of the model
 xmax                            = 4000.d0        ! abscissa of right side of the model
 nx                              = 80             ! number of elements along X
-nz                              = 60             ! number of elements along Z
 ngnod                           = 4              ! noeuds de controle pour blocs (4 ou 9)
-ratio                           = 0.967741       ! ratio pour separation en deux zones
 initialfield                    = .false.        ! use a plane wave as source or not
 ireadmodel                      = .false.        ! read external earth model or not
 TURN_ANISOTROPY_ON              = .false.        ! turn anisotropy on or off
@@ -39,7 +36,7 @@
 #
 # source parameters
 #
-isource_surf                    = .true.        ! sources dans le volume ou a la surface
+isource_surf                    = .true.         ! sources dans le volume ou a la surface
 xs                              = 2000.          ! source location x in meters
 zs                              = 1500.          ! source location z in meters
 f0                              = 10.0           ! central source frequency (Hz)
@@ -50,35 +47,36 @@
 #
 # receiver line parameters
 #
-ienreg_surf                     = .false.        ! enregistrement volume ou surface
-isismostype                     = 2              ! record 1=displ 2=veloc 3=accel
+ienreg_surf                     = .true.         ! enregistrement volume ou surface
+isismostype                     = 1              ! record 1=displ 2=veloc 3=accel
 nrec                            = 11             ! number of receivers
-xdeb                            = 0.             ! first receiver x in meters
+xdeb                            = 300.           ! first receiver x in meters
 zdeb                            = 2200.          ! first receiver z in meters
-xfin                            = 4000.          ! last receiver x in meters
+xfin                            = 3700.          ! last receiver x in meters
 zfin                            = 2200.          ! last receiver z in meters
 anglerec                        = 0.d0           ! angle to rotate components at receivers
 #
 # display parameters
 #
 itaff                           = 100            ! display frequency in time steps
-ivecttype                       = 2              ! display 1=displ 2=veloc 3=accel
+ivecttype                       = 1              ! display 1=displ 2=veloc 3=accel
 cutvect                         = 1.             ! amplitude min en % pour vector plots
 imeshvect                       = .true.         ! display mesh on vector plots or not
 imodelvect                      = .false.        ! display velocity model on vector plots
 iboundvect                      = .true.         ! display boundary conditions on plots
 interpol                        = .true.         ! interpolation of the display or not
 iptsdisp                        = 6              ! points for interpolation of display
-isubsamp                        = 2              ! subsampling of color snapshots
+isubsamp                        = 1              ! subsampling of color snapshots
 ignuplot                        = .false.        ! generate a GNUPLOT file for the grid
 ioutputgrid                     = .false.        ! save the grid in a text file or not
 #
 # velocity and density model (nx,nz)
 #
-nbmodels                        = 2              ! nb de modeles differents (rho,vp,vs)
-1 0 2200.d0 2500.d0 1443.375d0 0 0
-2 0 2200.d0 2500.d0 1443.375d0 0 0
-nbzone                          = 2              ! nb of zones and model number for each
-1 80 1 60 1
-3 51 3 38 2
-
+nbmodels                        = 3              ! nb de modeles differents (rho,vp,vs)
+1 0 2700.d0 3000.d0 1732.051d0 0 0
+2 0 2500.d0 2700.d0 1558.845d0 0 0
+3 0 2200.d0 2500.d0 1443.375d0 0 0
+nbzone                          = 3              ! nb of zones and model number for each
+1 80  1 20 1
+1 80 21 40 2
+1 80 41 60 3

Modified: seismo/2D/SPECFEM2D/trunk/MAILLE90/Par_cours_M2_UPPA
===================================================================
--- seismo/2D/SPECFEM2D/trunk/MAILLE90/Par_cours_M2_UPPA	2004-12-14 15:22:01 UTC (rev 8432)
+++ seismo/2D/SPECFEM2D/trunk/MAILLE90/Par_cours_M2_UPPA	2007-12-07 23:45:56 UTC (rev 8433)
@@ -6,20 +6,17 @@
 # ----------------------------------------------------------------
 <-                              ->
 #
-# File names and path for different outputs
+# title of job, and file that contains interface data
 #
 title                           = Test pour cours M2 UPPA
-topofile                        = topo_cours_M2_UPPA_curved.dat
-interffile                      = none
+interfacesfile                  = interfaces_cours_M2_UPPA_curved.dat
 #
 # geometry of the model (origin lower-left corner = 0,0) and mesh description
 #
-xmin                            = 0.0d0          ! abscissa of left side of the model
+xmin                            = 0.d0           ! abscissa of left side of the model
 xmax                            = 4000.d0        ! abscissa of right side of the model
 nx                              = 80             ! number of elements along X
-nz                              = 60             ! number of elements along Z
 ngnod                           = 4              ! noeuds de controle pour blocs (4 ou 9)
-ratio                           = 0.967741       ! ratio pour separation en deux zones
 initialfield                    = .false.        ! use a plane wave as source or not
 ireadmodel                      = .false.        ! read external earth model or not
 TURN_ANISOTROPY_ON              = .false.        ! turn anisotropy on or off
@@ -39,7 +36,7 @@
 #
 # source parameters
 #
-isource_surf                    = .true.        ! sources dans le volume ou a la surface
+isource_surf                    = .true.         ! sources dans le volume ou a la surface
 xs                              = 2000.          ! source location x in meters
 zs                              = 1500.          ! source location z in meters
 f0                              = 10.0           ! central source frequency (Hz)
@@ -50,35 +47,36 @@
 #
 # receiver line parameters
 #
-ienreg_surf                     = .false.        ! enregistrement volume ou surface
-isismostype                     = 2              ! record 1=displ 2=veloc 3=accel
+ienreg_surf                     = .true.         ! enregistrement volume ou surface
+isismostype                     = 1              ! record 1=displ 2=veloc 3=accel
 nrec                            = 11             ! number of receivers
-xdeb                            = 0.             ! first receiver x in meters
+xdeb                            = 300.           ! first receiver x in meters
 zdeb                            = 2200.          ! first receiver z in meters
-xfin                            = 4000.          ! last receiver x in meters
+xfin                            = 3700.          ! last receiver x in meters
 zfin                            = 2200.          ! last receiver z in meters
 anglerec                        = 0.d0           ! angle to rotate components at receivers
 #
 # display parameters
 #
 itaff                           = 100            ! display frequency in time steps
-ivecttype                       = 2              ! display 1=displ 2=veloc 3=accel
+ivecttype                       = 1              ! display 1=displ 2=veloc 3=accel
 cutvect                         = 1.             ! amplitude min en % pour vector plots
 imeshvect                       = .true.         ! display mesh on vector plots or not
 imodelvect                      = .false.        ! display velocity model on vector plots
 iboundvect                      = .true.         ! display boundary conditions on plots
 interpol                        = .true.         ! interpolation of the display or not
 iptsdisp                        = 6              ! points for interpolation of display
-isubsamp                        = 2              ! subsampling of color snapshots
+isubsamp                        = 1              ! subsampling of color snapshots
 ignuplot                        = .false.        ! generate a GNUPLOT file for the grid
 ioutputgrid                     = .false.        ! save the grid in a text file or not
 #
 # velocity and density model (nx,nz)
 #
-nbmodels                        = 2              ! nb de modeles differents (rho,vp,vs)
-1 0 2200.d0 2500.d0 1443.375d0 0 0
-2 0 2200.d0 2500.d0 1443.375d0 0 0
-nbzone                          = 2              ! nb of zones and model number for each
-1 80 1 60 1
-3 51 3 38 2
-
+nbmodels                        = 3              ! nb de modeles differents (rho,vp,vs)
+1 0 2700.d0 3000.d0 1732.051d0 0 0
+2 0 2500.d0 2700.d0 1558.845d0 0 0
+3 0 2200.d0 2500.d0 1443.375d0 0 0
+nbzone                          = 3              ! nb of zones and model number for each
+1 80  1 20 1
+1 80 21 40 2
+1 80 41 60 3

Deleted: seismo/2D/SPECFEM2D/trunk/MAILLE90/interf_paco.dat
===================================================================
--- seismo/2D/SPECFEM2D/trunk/MAILLE90/interf_paco.dat	2004-12-14 15:22:01 UTC (rev 8432)
+++ seismo/2D/SPECFEM2D/trunk/MAILLE90/interf_paco.dat	2007-12-07 23:45:56 UTC (rev 8433)
@@ -1,212 +0,0 @@
-   211
-   0.    4000.0000000000
-     33.333333333333    4000.0000000000
-     66.666666666667    4000.0000000000
-    100.000000000000    4000.0000000000
-     133.33333333333    4000.0000000000
-     166.66666666667    4000.0000000000
-     200.00000000000    4000.0000000000
-     233.33333333333    4000.0000000000
-     266.66666666667    4000.0000000000
-     300.00000000000    4000.0000000000
-     333.33333333333    4000.0000000000
-     366.66666666667    4000.0000000000
-     400.00000000000    4000.0000000000
-     433.33333333333    4000.0000000000
-     466.66666666667    4000.0000000000
-     500.00000000000    4000.0000000000
-     533.33333333333    4000.0000000000
-     566.66666666667    4000.0000000000
-     600.00000000000    4000.0000000000
-     633.33333333333    4000.0000000000
-     666.66666666667    4000.0000000000
-     700.00000000000    4000.0000000000
-     733.33333333333    4000.0000000000
-     766.66666666667    4000.0000000000
-     800.00000000000    4000.0000000000
-     833.33333333333    4000.0000000000
-     866.66666666667    4000.0000000000
-     900.00000000000    4000.0000000000
-     933.33333333333    4000.0000000000
-     966.66666666667    4000.0000000000
-    1000.00000000000    4000.0000000000
-     1033.3333333333    4000.0000000000
-     1066.6666666667    4000.0000000000
-     1100.0000000000    4000.0000000000
-     1133.3333333333    4000.0000000000
-     1166.6666666667    4000.0000000000
-     1200.0000000000    4000.0000000000
-     1233.3333333333    4000.0000000000
-     1266.6666666667    4000.0000000000
-     1300.0000000000    4000.0000000000
-     1333.3333333333    4000.0000000000
-     1366.6666666667    4000.0000000000
-     1400.0000000000    4000.0000000000
-     1433.3333333333    4000.0000000000
-     1466.6666666667    4000.0000000000
-     1500.0000000000    4000.0000000000
-     1533.3333333333    4000.0000000000
-     1566.6666666667    4000.0000000000
-     1600.0000000000    4000.0000000000
-     1633.3333333333    4000.0000000000
-     1666.6666666667    4000.0000000000
-     1700.0000000000    4000.0000000000
-     1733.3333333333    4000.0000000000
-     1766.6666666667    4000.0000000000
-     1800.0000000000    4000.0000000000
-     1833.3333333333    4000.0000000000
-     1866.6666666667    4000.0000000000
-     1900.0000000000    4000.0000000000
-     1933.3333333333    4000.0000000000
-     1966.6666666667    4000.0000000000
-     2000.0000000000    4000.0000000000
-     2033.3333333333    4000.0000000000
-     2066.6666666667    4000.0000000000
-     2100.0000000000    4000.0000000000
-     2133.3333333333    4000.0000000000
-     2166.6666666667    4000.0000000000
-     2200.0000000000    4000.0000000000
-     2233.3333333333    4000.0000000000
-     2266.6666666667    4000.0000000000
-     2300.0000000000    4000.0000000000
-     2333.3333333333    4000.0000000000
-     2366.6666666667    4000.0000000000
-     2400.0000000000    4000.0000000000
-     2433.3333333333    4000.0000000000
-     2466.6666666667    4000.0000000000
-     2500.0000000000    4000.0000000000
-     2533.3333333333    3998.6304738452
-     2566.6666666667    3994.5369001959
-     2600.0000000000    3987.7641291015
-     2633.3333333333    3978.3863644593
-     2666.6666666667    3966.5063510209
-     2700.0000000000    3952.2542486992
-     2733.3333333333    3935.7862065095
-     2766.6666666667    3917.2826517676
-     2800.0000000000    3896.9463132909
-     2833.3333333333    3875.0000002591
-     2866.6666666667    3851.6841610696
-     2900.0000000000    3827.2542489351
-     2933.3333333333    3801.9779230848
-     2966.6666666667    3776.1321162334
-     3000.0000000000    3750.0000004487
-     3033.3333333333    3723.8678846591
-     3066.6666666667    3698.0220777930
-     3100.0000000000    3672.7457519184
-     3133.3333333333    3648.3158397503
-     3166.6666666667    3625.0000005181
-     3200.0000000000    3603.0536874351
-     3233.3333333333    3582.7173488994
-     3266.6666666667    3564.2137940910
-     3300.0000000000    3547.7457518283
-     3333.3333333333    3533.4936494278
-     3366.6666666667    3521.6136359057
-     3400.0000000000    3512.2358711758
-     3433.3333333333    3505.4630999907
-     3466.6666666667    3501.3695262486
-     3500.0000000000    3500.0000000000
-     3533.3333333333    3501.3695260610
-     3566.6666666667    3505.4630996175
-     3600.0000000000    3512.2358706212
-     3633.3333333333    3521.6136351757
-     3666.6666666667    3533.4936485304
-     3700.0000000000    3547.7457507733
-     3733.3333333333    3564.2137928900
-     3766.6666666667    3582.7173475655
-     3800.0000000000    3603.0536859830
-     3833.3333333333    3624.9999989637
-     3866.6666666667    3648.3158381106
-     3900.0000000000    3672.7457502113
-     3933.3333333333    3698.0220760373
-     3966.6666666667    3723.8678828740
-     4000.0000000000    3749.9999986538
-     4033.3333333333    3776.1321144484
-     4066.6666666667    3801.9779213292
-     4100.0000000000    3827.2542472281
-     4133.3333333333    3851.6841594298
-     4166.6666666667    3874.9999987046
-     4200.0000000000    3896.9463118388
-     4233.3333333333    3917.2826504337
-     4266.6666666667    3935.7862053084
-     4300.0000000000    3952.2542476442
-     4333.3333333333    3966.5063501234
-     4366.6666666667    3978.3863637293
-     4400.0000000000    3987.7641285469
-     4433.3333333333    3994.5368998227
-     4466.6666666667    3998.6304736576
-     4500.0000000000    4000.0000000000
-     4533.3333333333    4000.0000000000
-     4566.6666666667    4000.0000000000
-     4600.0000000000    4000.0000000000
-     4633.3333333333    4000.0000000000
-     4666.6666666667    4000.0000000000
-     4700.0000000000    4000.0000000000
-     4733.3333333333    4000.0000000000
-     4766.6666666667    4000.0000000000
-     4800.0000000000    4000.0000000000
-     4833.3333333333    4000.0000000000
-     4866.6666666667    4000.0000000000
-     4900.0000000000    4000.0000000000
-     4933.3333333333    4000.0000000000
-     4966.6666666667    4000.0000000000
-     5000.0000000000    4000.0000000000
-     5033.3333333333    4000.0000000000
-     5066.6666666667    4000.0000000000
-     5100.0000000000    4000.0000000000
-     5133.3333333333    4000.0000000000
-     5166.6666666667    4000.0000000000
-     5200.0000000000    4000.0000000000
-     5233.3333333333    4000.0000000000
-     5266.6666666667    4000.0000000000
-     5300.0000000000    4000.0000000000
-     5333.3333333333    4000.0000000000
-     5366.6666666667    4000.0000000000
-     5400.0000000000    4000.0000000000
-     5433.3333333333    4000.0000000000
-     5466.6666666667    4000.0000000000
-     5500.0000000000    4000.0000000000
-     5533.3333333333    4000.0000000000
-     5566.6666666667    4000.0000000000
-     5600.0000000000    4000.0000000000
-     5633.3333333333    4000.0000000000
-     5666.6666666667    4000.0000000000
-     5700.0000000000    4000.0000000000
-     5733.3333333333    4000.0000000000
-     5766.6666666667    4000.0000000000
-     5800.0000000000    4000.0000000000
-     5833.3333333333    4000.0000000000
-     5866.6666666667    4000.0000000000
-     5900.0000000000    4000.0000000000
-     5933.3333333333    4000.0000000000
-     5966.6666666667    4000.0000000000
-     6000.0000000000    4000.0000000000
-     6033.3333333333    4000.0000000000
-     6066.6666666667    4000.0000000000
-     6100.0000000000    4000.0000000000
-     6133.3333333333    4000.0000000000
-     6166.6666666667    4000.0000000000
-     6200.0000000000    4000.0000000000
-     6233.3333333333    4000.0000000000
-     6266.6666666667    4000.0000000000
-     6300.0000000000    4000.0000000000
-     6333.3333333333    4000.0000000000
-     6366.6666666667    4000.0000000000
-     6400.0000000000    4000.0000000000
-     6433.3333333333    4000.0000000000
-     6466.6666666667    4000.0000000000
-     6500.0000000000    4000.0000000000
-     6533.3333333333    4000.0000000000
-     6566.6666666667    4000.0000000000
-     6600.0000000000    4000.0000000000
-     6633.3333333333    4000.0000000000
-     6666.6666666667    4000.0000000000
-     6700.0000000000    4000.0000000000
-     6733.3333333333    4000.0000000000
-     6766.6666666667    4000.0000000000
-     6800.0000000000    4000.0000000000
-     6833.3333333333    4000.0000000000
-     6866.6666666667    4000.0000000000
-     6900.0000000000    4000.0000000000
-     6933.3333333333    4000.0000000000
-     6966.6666666667    4000.0000000000
-     7000.0000000000    4000.0000000000

Modified: seismo/2D/SPECFEM2D/trunk/MAILLE90/meshfem2D.f90
===================================================================
--- seismo/2D/SPECFEM2D/trunk/MAILLE90/meshfem2D.f90	2004-12-14 15:22:01 UTC (rev 8432)
+++ seismo/2D/SPECFEM2D/trunk/MAILLE90/meshfem2D.f90	2007-12-07 23:45:56 UTC (rev 8433)
@@ -1,19 +1,19 @@
 
 !========================================================================
 !
-!                   M E S H F E M 2 D  Version 5.0
+!                   M E S H F E M 2 D  Version 5.1
 !                   ------------------------------
 !
 !                         Dimitri Komatitsch
 !          Universite de Pau et des Pays de l'Adour, France
 !
-!                          (c) May 2004
+!                         (c) December 2004
 !
 !========================================================================
 
 !========================================================================
 !
-!  Mesh generator for SPECFEM2D version 5.0
+!  Basic mesh generator for SPECFEM2D
 !
 !========================================================================
 
@@ -23,38 +23,42 @@
 
 ! definir les tableaux pour allocation dynamique
 
-! coordinates of the grid points
-  double precision, allocatable :: x(:,:),z(:,:)
+  integer, parameter :: ANISOTROPIC_MATERIAL = 1
 
-! variables needed to compute the transformation
-  double precision, allocatable :: psi(:),eta(:),absx(:), &
-      a00(:),a01(:),valeta(:),bot0(:),top0(:)
+! coordinates of the grid points of the mesh
+  double precision, dimension(:,:), allocatable :: x,z
 
+! to compute the coordinate transformation
+  integer :: ioffset
+  double precision :: eta,absx,a00,a01,bot0,top0
+
 ! stockage du modele de vitesse et densite
-  double precision, allocatable :: rho(:),cp(:),cs(:),aniso3(:),aniso4(:)
-  integer, allocatable :: icodemat(:)
-  integer, allocatable :: num_modele(:,:)
+  double precision, dimension(:), allocatable :: rho,cp,cs,aniso3,aniso4
+  integer, dimension(:), allocatable :: icodemat
+  integer, dimension(:,:), allocatable :: num_modele
 
-! the topography data
-  double precision, allocatable :: xtopo(:),ztopo(:),coefs_topo(:), &
-      xinterf(:),zinterf(:),coefs_interf(:)
+! interface data
+  integer interface_current,ipoint_current,number_of_interfaces,npoints_interface_bottom,npoints_interface_top
+  integer ilayer,number_of_layers,max_npoints_interface
+  double precision xinterface_dummy,zinterface_dummy,xinterface_dummy_previous
+  integer, dimension(:), allocatable :: nz_layer
+  double precision, dimension(:), allocatable :: &
+         xinterface_bottom,zinterface_bottom,coefs_interface_bottom, &
+         xinterface_top,zinterface_top,coefs_interface_top
 
 ! for the source
-  double precision xs,zs,f0,t0,angle,factor
   integer isource_type
+  double precision xs,zs,f0,t0,angle,factor
 
 ! arrays for the receivers
-  double precision, allocatable :: xrec(:),zrec(:)
+  double precision, dimension(:), allocatable :: xrec,zrec
 
-! nom du fichier GNUPLOT contenant la grille
-  character(len=50) file1
-
-  character(len=50) interffile,topofile,title
+  character(len=50) interfacesfile,title
   character(len=34) junk
 
-  integer imatnum,inumabs,inumelem,netyp
-  integer nelemabs,npgeo,nspec,ninterf,ntopo
-  integer k,icol,ili,istepx,istepz,ncut,ix,iz,irec,i,j
+  integer imatnum,inumabs,inumelem
+  integer nelemabs,npgeo,nspec
+  integer k,icol,ili,istepx,istepz,ix,iz,irec,i,j
   integer ixdebzone,ixfinzone,izdebzone,izfinzone,imodnum
   integer izone,imodele,nbzone,nbmodeles
   integer itaff,iptsdisp,isubsamp,nrec
@@ -64,30 +68,24 @@
 
   logical codehaut,codebas,codegauche,codedroite
 
-  double precision ratio
   double precision tang1,tangN,vpzone,vszone
-  double precision cutvect
-  double precision xspacerec,zspacerec
-  double precision anglerec,xfin,zfin,xdeb,zdeb,xmin,xmax
-  double precision dt
+  double precision cutvect,xspacerec,zspacerec
+  double precision anglerec,xfin,zfin,xdeb,zdeb,xmin,xmax,dt
   double precision rhoread,cpread,csread,aniso3read,aniso4read
 
   logical interpol,ignuplot,ireadmodel,ioutputgrid
   logical abshaut,absbas,absgauche,absdroite
-  logical isource_surf,ienreg_surf
-  logical imeshvect
-  logical initialfield
-  logical imodelvect,iboundvect
+  logical isource_surf,ienreg_surf,imeshvect,initialfield,imodelvect,iboundvect
   logical TURN_ANISOTROPY_ON,TURN_ATTENUATION_ON
 
   integer, external :: num
-  double precision, external :: bottom,spl,dens
+  double precision, external :: value_spline
 
 ! ***
 ! *** read the parameter file
 ! ***
 
-  print *,' Reading the parameter file ... '
+  print *,'Reading the parameter file ... '
   print *
 
   open(unit=10,file='Par',status='old')
@@ -104,9 +102,8 @@
   enddo
 
 ! read file names and path for output
-  read(10,3)junk,title
-  read(10,3)junk,topofile
-  read(10,3)junk,interffile
+  read(10,3) junk,title
+  read(10,3) junk,interfacesfile
 
   write(*,*) 'Titre de la simulation'
   write(*,*) title
@@ -118,17 +115,92 @@
   read(10,*)
 
 ! read grid parameters
-  read(10,1)junk,xmin
-  read(10,1)junk,xmax
-  read(10,2)junk,nx
-  read(10,2)junk,nz
-  read(10,2)junk,ngnod
-  read(10,1)junk,ratio
-  read(10,4)junk,initialfield
-  read(10,4)junk,ireadmodel
-  read(10,4)junk,TURN_ANISOTROPY_ON
-  read(10,4)junk,TURN_ATTENUATION_ON
+  read(10,1) junk,xmin
+  read(10,1) junk,xmax
+  read(10,2) junk,nx
+  read(10,2) junk,ngnod
+  read(10,4) junk,initialfield
+  read(10,4) junk,ireadmodel
+  read(10,4) junk,TURN_ANISOTROPY_ON
+  read(10,4) junk,TURN_ATTENUATION_ON
 
+! get interface data from external file to count the spectral elements along Z
+  print *,'Reading interface data from file ',interfacesfile(1:len_trim(interfacesfile)),' to count the spectral elements'
+  open(unit=15,file=interfacesfile,status='old')
+
+  max_npoints_interface = -1
+
+! skip header
+  read(15,*)
+  read(15,*)
+  read(15,*)
+
+! read number of interfaces
+  read(15,*) number_of_interfaces
+  if(number_of_interfaces < 2) stop 'not enough interfaces (minimum is 2)'
+
+! skip header
+  read(15,*)
+  read(15,*)
+  read(15,*)
+
+! loop on all the interfaces
+  do interface_current = 1,number_of_interfaces
+
+! skip header
+    read(15,*)
+    read(15,*)
+    read(15,*)
+
+    read(15,*) npoints_interface_bottom
+    if(npoints_interface_bottom < 2) stop 'not enough interface points (minimum is 2)'
+    max_npoints_interface = max(npoints_interface_bottom,max_npoints_interface)
+    print *,'Reading ',npoints_interface_bottom,' points for interface ',interface_current
+
+! loop on all the points describing this interface
+    do ipoint_current = 1,npoints_interface_bottom
+      read(15,*) xinterface_dummy,zinterface_dummy
+      if(ipoint_current > 1 .and. xinterface_dummy <= xinterface_dummy_previous) &
+        stop 'interface points must be sorted in increasing X'
+      xinterface_dummy_previous = xinterface_dummy
+    enddo
+
+  enddo
+
+! define number of layers
+  number_of_layers = number_of_interfaces - 1
+
+  allocate(nz_layer(number_of_layers))
+
+! skip header
+    read(15,*)
+    read(15,*)
+    read(15,*)
+
+! loop on all the layers
+  do ilayer = 1,number_of_layers
+
+! skip header
+    read(15,*)
+    read(15,*)
+    read(15,*)
+
+! read number of spectral elements in vertical direction in this layer
+    read(15,*) nz_layer(ilayer)
+    if(nz_layer(ilayer) < 1) stop 'not enough spectral elements along Z in layer (minimum is 1)'
+    print *,'There are ',nz_layer(ilayer),' spectral elements along Z in layer ',ilayer
+
+  enddo
+
+  close(15)
+
+! compute total number of spectral elements in vertical direction
+  nz = sum(nz_layer)
+
+  print *
+  print *,'Total number of spectral elements along Z = ',nz
+  print *
+
   nxread = nx
   nzread = nz
 
@@ -144,10 +216,10 @@
   read(10,*)
 
 ! read absorbing boundaries parameters
-  read(10,4)junk,abshaut
-  read(10,4)junk,absbas
-  read(10,4)junk,absgauche
-  read(10,4)junk,absdroite
+  read(10,4) junk,abshaut
+  read(10,4) junk,absbas
+  read(10,4) junk,absgauche
+  read(10,4) junk,absdroite
 
 ! skip comment
   read(10,*)
@@ -155,8 +227,8 @@
   read(10,*)
 
 ! read time step parameters
-  read(10,2)junk,nt
-  read(10,1)junk,dt
+  read(10,2) junk,nt
+  read(10,1) junk,dt
 
 ! skip comment
   read(10,*)
@@ -164,20 +236,20 @@
   read(10,*)
 
 ! read source parameters
-  read(10,4)junk,isource_surf
-  read(10,1)junk,xs
-  read(10,1)junk,zs
-  read(10,1)junk,f0
-  read(10,1)junk,t0
-  read(10,2)junk,isource_type
-  read(10,1)junk,angle
-  read(10,1)junk,factor
+  read(10,4) junk,isource_surf
+  read(10,1) junk,xs
+  read(10,1) junk,zs
+  read(10,1) junk,f0
+  read(10,1) junk,t0
+  read(10,2) junk,isource_type
+  read(10,1) junk,angle
+  read(10,1) junk,factor
 
   print *
   print *,'Source:'
   print *,'Position xs, zs = ',xs,zs
   print *,'Frequency, delay = ',f0,t0
-  print *,'Source type (1=force 2=explo) : ',isource_type
+  print *,'Source type (1=force 2=explosion) : ',isource_type
   print *,'Angle of the source if force = ',angle
   print *,'Multiplying factor = ',factor
 
@@ -187,22 +259,22 @@
   read(10,*)
 
 ! read receivers line parameters
-  read(10,4)junk,ienreg_surf
-  read(10,2)junk,isismostype
-  read(10,2)junk,nrec
-  read(10,1)junk,xdeb
-  read(10,1)junk,zdeb
-  read(10,1)junk,xfin
-  read(10,1)junk,zfin
-  read(10,1)junk,anglerec
+  read(10,4) junk,ienreg_surf
+  read(10,2) junk,isismostype
+  read(10,2) junk,nrec
+  read(10,1) junk,xdeb
+  read(10,1) junk,zdeb
+  read(10,1) junk,xfin
+  read(10,1) junk,zfin
+  read(10,1) junk,anglerec
 
   allocate(xrec(nrec))
   allocate(zrec(nrec))
 
   print *
   print *,'There are ',nrec,' receivers'
-  xspacerec=(xfin-xdeb)/dble(nrec-1)
-  zspacerec=(zfin-zdeb)/dble(nrec-1)
+  xspacerec = (xfin-xdeb) / dble(nrec-1)
+  zspacerec = (zfin-zdeb) / dble(nrec-1)
   do i=1,nrec
     xrec(i) = xdeb + dble(i-1)*xspacerec
     zrec(i) = zdeb + dble(i-1)*zspacerec
@@ -214,17 +286,17 @@
   read(10,*)
 
 ! read display parameters
-  read(10,2)junk,itaff
-  read(10,2)junk,ivecttype
-  read(10,1)junk,cutvect
-  read(10,4)junk,imeshvect
-  read(10,4)junk,imodelvect
-  read(10,4)junk,iboundvect
-  read(10,4)junk,interpol
-  read(10,2)junk,iptsdisp
-  read(10,2)junk,isubsamp
-  read(10,4)junk,ignuplot
-  read(10,4)junk,ioutputgrid
+  read(10,2) junk,itaff
+  read(10,2) junk,ivecttype
+  read(10,1) junk,cutvect
+  read(10,4) junk,imeshvect
+  read(10,4) junk,imodelvect
+  read(10,4) junk,iboundvect
+  read(10,4) junk,interpol
+  read(10,2) junk,iptsdisp
+  read(10,2) junk,isubsamp
+  read(10,4) junk,ignuplot
+  read(10,4) junk,ioutputgrid
 
 ! skip comment
   read(10,*)
@@ -232,8 +304,7 @@
   read(10,*)
 
 ! lecture des differents modeles de materiaux
-
-  read(10,2)junk,nbmodeles
+  read(10,2) junk,nbmodeles
   if(nbmodeles <= 0) stop 'Negative number of models not allowed !!'
 
   allocate(icodemat(nbmodeles))
@@ -253,24 +324,23 @@
   num_modele(:,:) = 0
 
   do imodele=1,nbmodeles
-      read(10,*) i,icodematread,rhoread,cpread,csread,aniso3read,aniso4read
-      if(i<1 .or. i>nbmodeles) stop 'Wrong material point number'
-      icodemat(i) = icodematread
-      rho(i) = rhoread
-      cp(i) = cpread
-      cs(i) = csread
-      aniso3(i) = aniso3read
-      aniso4(i) = aniso4read
-      if(i <= 0) stop 'Negative model number not allowed !!'
-      if (rho(i) < 0.d0 .or. cp(i) < 0.d0 .or. cs(i) < 0.d0) &
-          stop 'Negative value of velocity or density'
+    read(10,*) i,icodematread,rhoread,cpread,csread,aniso3read,aniso4read
+    if(i < 1 .or. i > nbmodeles) stop 'Wrong material point number'
+    icodemat(i) = icodematread
+    rho(i) = rhoread
+    cp(i) = cpread
+    cs(i) = csread
+    aniso3(i) = aniso3read
+    aniso4(i) = aniso4read
+    if(i <= 0) stop 'Negative model number not allowed !!'
+    if(rho(i) <= 0.d0 .or. cp(i) <= 0.d0 .or. cs(i) < 0.d0) stop 'negative value of velocity or density'
   enddo
 
   print *
   print *, 'Nb de modeles de roche = ',nbmodeles
   print *
   do i=1,nbmodeles
-    if(icodemat(i) /= 2) then
+    if(icodemat(i) /= ANISOTROPIC_MATERIAL) then
       print *,'Modele #',i,' isotrope'
       print *,'rho,cp,cs = ',rho(i),cp(i),cs(i)
     else
@@ -280,8 +350,8 @@
   enddo
   print *
 
-! *** lecture des numeros de modele des differentes zones
-  read(10,2)junk,nbzone
+! lecture des numeros de modele des differentes zones
+  read(10,2) junk,nbzone
 
   if(nbzone <= 0) stop 'Negative number of zones not allowed !!'
 
@@ -289,54 +359,56 @@
   print *, 'Nb de zones du modele = ',nbzone
   print *
 
-  do izone=1,nbzone
-      read(10,*) ixdebzone,ixfinzone,izdebzone,izfinzone,imodnum
+  do izone = 1,nbzone
 
-  if (imodnum < 1) stop 'Negative model number not allowed !!'
-  if (ixdebzone < 1) stop 'Left coordinate of zone negative !!'
-  if (ixfinzone > nxread) stop 'Right coordinate of zone too high !!'
-  if (izdebzone < 1) stop 'Bottom coordinate of zone negative !!'
-  if (izfinzone > nzread) stop 'Top coordinate of zone too high !!'
+    read(10,*) ixdebzone,ixfinzone,izdebzone,izfinzone,imodnum
 
-      print *,'Zone ',izone
-      print *,'IX from ',ixdebzone,' to ',ixfinzone
-      print *,'IZ from ',izdebzone,' to ',izfinzone
-  if(icodemat(imodnum) /= 2) then
-      vpzone = cp(imodnum)
-      vszone = cs(imodnum)
-      print *,'Model # ',imodnum,' isotrope'
-      print *,'vp = ',vpzone
-      print *,'vs = ',vszone
-      print *,'rho = ',rho(imodnum)
-      print *,'Poisson''s ratio = ', &
-        0.5d0*(vpzone*vpzone-2.d0*vszone*vszone)/(vpzone*vpzone-vszone*vszone)
+    if(imodnum < 1) stop 'Negative model number not allowed !!'
+    if(ixdebzone < 1) stop 'Left coordinate of zone negative !!'
+    if(ixfinzone > nxread) stop 'Right coordinate of zone too high !!'
+    if(izdebzone < 1) stop 'Bottom coordinate of zone negative !!'
+    if(izfinzone > nzread) stop 'Top coordinate of zone too high !!'
+
+    print *,'Zone ',izone
+    print *,'IX from ',ixdebzone,' to ',ixfinzone
+    print *,'IZ from ',izdebzone,' to ',izfinzone
+
+  if(icodemat(imodnum) /= ANISOTROPIC_MATERIAL) then
+    vpzone = cp(imodnum)
+    vszone = cs(imodnum)
+    print *,'Model # ',imodnum,' isotrope'
+    print *,'vp = ',vpzone
+    print *,'vs = ',vszone
+    print *,'rho = ',rho(imodnum)
+    print *,'Poisson''s ratio = ', &
+      0.5d0*(vpzone*vpzone-2.d0*vszone*vszone) / (vpzone*vpzone-vszone*vszone)
   else
-      print *,'Model # ',imodnum,' anisotrope'
-      print *,'c11 = ',cp(imodnum)
-      print *,'c13 = ',cs(imodnum)
-      print *,'c33 = ',aniso3(imodnum)
-      print *,'c44 = ',aniso4(imodnum)
-      print *,'rho = ',rho(imodnum)
+    print *,'Model # ',imodnum,' anisotrope'
+    print *,'c11 = ',cp(imodnum)
+    print *,'c13 = ',cs(imodnum)
+    print *,'c33 = ',aniso3(imodnum)
+    print *,'c44 = ',aniso4(imodnum)
+    print *,'rho = ',rho(imodnum)
   endif
-      print *,' -----'
+  print *,' -----'
 
 ! stocker le modele de vitesse et densite
-   do i=ixdebzone,ixfinzone
-      do j=izdebzone,izfinzone
-        if(ngnod == 4) then
-          num_modele(i,j) = imodnum
-        else
-          num_modele(2*(i-1)+1,2*(j-1)+1) = imodnum
-          num_modele(2*(i-1)+1,2*(j-1)+2) = imodnum
-          num_modele(2*(i-1)+2,2*(j-1)+1) = imodnum
-          num_modele(2*(i-1)+2,2*(j-1)+2) = imodnum
-        endif
-      enddo
+   do i = ixdebzone,ixfinzone
+     do j = izdebzone,izfinzone
+       if(ngnod == 4) then
+         num_modele(i,j) = imodnum
+       else
+         num_modele(2*(i-1)+1,2*(j-1)+1) = imodnum
+         num_modele(2*(i-1)+1,2*(j-1)+2) = imodnum
+         num_modele(2*(i-1)+2,2*(j-1)+1) = imodnum
+         num_modele(2*(i-1)+2,2*(j-1)+2) = imodnum
+       endif
+     enddo
    enddo
 
   enddo
 
-  if (minval(num_modele) <= 0) stop 'Velocity model not entirely set...'
+  if(minval(num_modele) <= 0) stop 'Velocity model not entirely set...'
 
   close(10)
 
@@ -345,28 +417,8 @@
 
 ! --------- fin lecture fichier parametres --------------
 
-  allocate(psi(0:nx))
-  allocate(eta(0:nz))
-  allocate(absx(0:nx))
-  allocate(a00(0:nz))
-  allocate(a01(0:nz))
-  allocate(valeta(0:nz))
-  allocate(bot0(0:nx))
-  allocate(top0(0:nx))
+  if(ngnod /= 4 .and. ngnod /= 9) stop 'erreur ngnod different de 4 ou 9 !!'
 
-! calcul des points regulierement espaces
-  do i=0,nx
-    psi(i) = i/dble(nx)
-  enddo
-
-  do j=0,nz
-    eta(j) = j/dble(nz)
-  enddo
-
-! quelques verifications de base a faire
-
-  if(ngnod /= 4.and.ngnod /= 9) stop 'erreur ngnod different de 4 ou 9 !!'
-
   print *
   if(ngnod == 4) then
     print *,'Le maillage comporte ',nx,' x ',nz,' elements'
@@ -377,211 +429,167 @@
   print *,'Les elements de controle sont des elements ',ngnod,' noeuds'
   print *
 
-!------------------------------------------------------
+!---
 
+! allocate arrays for the grid
   allocate(x(0:nx,0:nz))
   allocate(z(0:nx,0:nz))
 
-  x(:,:)=0.d0
-  z(:,:)=0.d0
+  x(:,:) = 0.d0
+  z(:,:) = 0.d0
 
-! get topography data from external file
-  print *,'Reading topography from file ',topofile
-  open(unit=15,file=topofile,status='old')
-  read(15,*) ntopo
-  if (ntopo < 2) stop 'Not enough topography points (min 2)'
-  print *,'Reading ',ntopo,' points from topography file'
-  print *
+! get interface data from external file
+  print *,'Reading interface data from file ',interfacesfile(1:len_trim(interfacesfile))
+  open(unit=15,file=interfacesfile,status='old')
 
-  allocate(xtopo(ntopo))
-  allocate(ztopo(ntopo))
-  allocate(coefs_topo(ntopo))
+  allocate(xinterface_bottom(max_npoints_interface))
+  allocate(zinterface_bottom(max_npoints_interface))
+  allocate(coefs_interface_bottom(max_npoints_interface))
 
-  do i=1,ntopo
-    read(15,*) xtopo(i),ztopo(i)
-  enddo
-  close(15)
+  allocate(xinterface_top(max_npoints_interface))
+  allocate(zinterface_top(max_npoints_interface))
+  allocate(coefs_interface_top(max_npoints_interface))
 
-! get interface data from external file, if any
-  if(interffile /= 'none') then
-  print *,'Reading interface from file ',interffile
-  open(unit=15,file=interffile,status='old')
-  read(15,*) ninterf
-  if (ninterf < 2) stop 'Not enough interface points (min 2)'
-  print *,'Reading ',ninterf,' points from interface file'
+! skip header
+  read(15,*)
+  read(15,*)
+  read(15,*)
 
-  allocate(xinterf(ninterf))
-  allocate(zinterf(ninterf))
-  allocate(coefs_interf(ninterf))
+! read number of interfaces
+  read(15,*) number_of_interfaces
 
-  do i=1,ninterf
-    read(15,*) xinterf(i),zinterf(i)
-  enddo
+! skip header
+  read(15,*)
+  read(15,*)
+  read(15,*)
 
-  close(15)
-  else
-    print *,'*** No interface file specified ***'
-  endif
+! skip header
+  read(15,*)
+  read(15,*)
+  read(15,*)
 
-! check the values read
-  print *
-  print *, 'Topography data points (x,z)'
-  print *, '----------------------------'
-  print *
-  print *, 'Topo 1 = (',xtopo(1),',',ztopo(1),')'
-  print *, 'Topo ntopo = (',xtopo(ntopo),',',ztopo(ntopo),')'
+! read bottom interface
+  read(15,*) npoints_interface_bottom
 
-!--- calculate the spline function for the topography
-!--- imposer les tangentes aux deux bords
-  tang1 = (ztopo(2)-ztopo(1))/(xtopo(2)-xtopo(1))
-  tangN = (ztopo(ntopo)-ztopo(ntopo-1))/(xtopo(ntopo)-xtopo(ntopo-1))
-  call spline(xtopo,ztopo,ntopo,tang1,tangN,coefs_topo)
+! loop on all the points describing this interface
+  do ipoint_current = 1,npoints_interface_bottom
+    read(15,*) xinterface_bottom(ipoint_current),zinterface_bottom(ipoint_current)
+  enddo
 
-!--- calculate the spline function for the interface
-!--- imposer les tangentes aux deux bords
-  if (interffile /= 'none') then
-  tang1 = (zinterf(2)-zinterf(1))/(xinterf(2)-xinterf(1))
-  tangN = (zinterf(ntopo)-zinterf(ntopo-1))/(xinterf(ntopo)-xinterf(ntopo-1))
-  call spline(xinterf,zinterf,Ninterf,tang1,tangN,coefs_interf)
-  endif
+! boucle sur toutes les couches
+  do ilayer = 1,number_of_layers
 
-! *** afficher limites du modele lu
-  print *
-  print *, 'Limites absolues modele fichier topo :'
-  print *
-  print *, 'Xmin = ',minval(xtopo),'   Xmax = ',maxval(xtopo)
-  print *, 'Zmin = ',minval(ztopo),'   Zmax = ',maxval(ztopo)
-  print *
+! skip header
+  read(15,*)
+  read(15,*)
+  read(15,*)
 
-! *** eventuellement modifier sources si sources en surface
-  print *
-  print *, 'Position (x,z) de la source'
-  print *
-  if(isource_surf) zs = spl(xs,xtopo,ztopo,coefs_topo,ntopo)
-  print *, 'Source = ',xs,zs
+! read top interface
+  read(15,*) npoints_interface_top
 
-! *** eventuellement modifier recepteurs si enregistrement en surface
-  print *
-  print *, 'Position (x,z) des ',nrec,' receivers'
-  print *
-  do irec=1,nrec
-   if(ienreg_surf) zrec(irec) = spl(xrec(irec),xtopo,ztopo,coefs_topo,ntopo)
-   print *, 'Receiver ',irec,' = ',xrec(irec),zrec(irec)
+! loop on all the points describing this interface
+  do ipoint_current = 1,npoints_interface_top
+    read(15,*) xinterface_top(ipoint_current),zinterface_top(ipoint_current)
   enddo
 
-!--- definition du maillage suivant X
-  do ix=0,nx
-    absx(ix) = dens(ix,psi,xmin,xmax,nx)
-  enddo
+! calculer le spline pour l'interface du bas, imposer la tangente aux deux bords
+  tang1 = (zinterface_bottom(2)-zinterface_bottom(1)) / (xinterface_bottom(2)-xinterface_bottom(1))
+  tangN = (zinterface_bottom(npoints_interface_bottom)-zinterface_bottom(npoints_interface_bottom-1)) / &
+          (xinterface_bottom(npoints_interface_bottom)-xinterface_bottom(npoints_interface_bottom-1))
+  call spline(xinterface_bottom,zinterface_bottom,npoints_interface_bottom,tang1,tangN,coefs_interface_bottom)
 
-  if (interffile == 'none') then
+! calculer le spline pour l'interface du haut, imposer la tangente aux deux bords
+  tang1 = (zinterface_top(2)-zinterface_top(1)) / (xinterface_top(2)-xinterface_top(1))
+  tangN = (zinterface_top(npoints_interface_top)-zinterface_top(npoints_interface_top-1)) / &
+          (xinterface_top(npoints_interface_top)-xinterface_top(npoints_interface_top-1))
+  call spline(xinterface_top,zinterface_top,npoints_interface_top,tang1,tangN,coefs_interface_top)
 
-! *** une seule zone si pas d'interface specifiee
+! tester si on est sur la derniere couche, qui contient la topographie
+  if(ilayer == number_of_layers) then
 
-  do iz=0,nz
-    valeta(iz) = eta(iz)
-    a00(iz) = 1-valeta(iz)
-    a01(iz) = valeta(iz)
-  enddo
+! modifier position de la source si source exactement en surface
+    if(isource_surf) zs = value_spline(xs,xinterface_top,zinterface_top,coefs_interface_top,npoints_interface_top)
 
-  do ix=0,nx
-    bot0(ix) = bottom(absx(ix))
-    top0(ix) = spl(absx(ix),xtopo,ztopo,coefs_topo,ntopo)
-  enddo
+! modifier position des recepteurs si enregistrement exactement en surface
+    if(ienreg_surf) then
+      do irec=1,nrec
+        zrec(irec) = value_spline(xrec(irec),xinterface_top,zinterface_top,coefs_interface_top,npoints_interface_top)
+      enddo
+    endif
 
-! valeurs de x et y pour display domaine physique
-  do ix=0,nx
-    do iz=0,nz
-      x(ix,iz) = absx(ix)
-      z(ix,iz) = a00(iz)*bot0(ix) + a01(iz)*top0(ix)
-    enddo
-  enddo
+  endif
 
+! calcul de l'offset de cette couche en nombre d'elements spectraux suivant Z
+  if(ilayer > 1) then
+    ioffset = sum(nz_layer(1:ilayer-1))
   else
+    ioffset = 0
+  endif
 
-! *** deux zones si topo
+!--- definition du maillage
 
-  ncut = nint(nz*ratio)
+    do ix = 0,nx
 
-! *** ZONE DU BAS ***
+! points regulierement espaces suivant X
+      absx = xmin + (xmax - xmin) * dble(ix) / dble(nx)
 
-  do iz=0,ncut
-    valeta(iz) = dble(iz)/(nz*ratio)
-    a00(iz) = 1-valeta(iz)
-    a01(iz) = valeta(iz)
-  enddo
+! value of the bottom and top splines
+      bot0 = value_spline(absx,xinterface_bottom,zinterface_bottom,coefs_interface_bottom,npoints_interface_bottom)
+      top0 = value_spline(absx,xinterface_top,zinterface_top,coefs_interface_top,npoints_interface_top)
 
-  do ix=0,nx
-    bot0(ix) = bottom(absx(ix))
-    top0(ix) = spl(absx(ix),xinterf,zinterf,coefs_interf,ninterf)
-  enddo
+      do iz = 0,nz_layer(ilayer)
 
-! valeurs de x et y pour display domaine physique
-  do ix=0,nx
-    do iz=0,ncut
-      x(ix,iz) = absx(ix)
-      z(ix,iz) = a00(iz)*bot0(ix) + a01(iz)*top0(ix)
-    enddo
-  enddo
+! linear interpolation between bottom and top
+        eta = dble(iz) / dble(nz_layer(ilayer))
+        a00 = 1.d0 - eta
+        a01 = eta
 
-! *** ZONE DU HAUT ***
+! coordinates of the grid points
+        x(ix,iz + ioffset) = absx
+        z(ix,iz + ioffset) = a00*bot0 + a01*top0
 
-  do iz=0,nz-ncut
-    valeta(iz) = dble(iz)/(nz*(1.d0-ratio))
-    a00(iz) = 1-valeta(iz)
-    a01(iz) = valeta(iz)
-  enddo
+      enddo
 
-  do ix=0,nx
-    bot0(ix) = spl(absx(ix),xinterf,zinterf,coefs_interf,ninterf)
-    top0(ix) = spl(absx(ix),xtopo,ztopo,coefs_topo,ntopo)
-  enddo
-
-! valeurs de x et y pour display domaine physique
-  do ix=0,nx
-    do iz=0,nz-ncut
-      x(ix,iz+ncut) = absx(ix)
-      z(ix,iz+ncut) = a00(iz)*bot0(ix) + a01(iz)*top0(ix)
     enddo
+
+! l'interface du haut devient celle du bas pour passer a la couche suivante
+    npoints_interface_bottom = npoints_interface_top
+    xinterface_bottom(:) = xinterface_top(:)
+    zinterface_bottom(:) = zinterface_top(:)
+
   enddo
 
-  endif
+  close(15)
 
 ! calculer min et max de X et Z sur la grille
   print *
-  print *, 'Valeurs min et max de X sur le maillage = ',minval(x),maxval(x)
-  print *, 'Valeurs min et max de Z sur le maillage = ',minval(z),maxval(z)
+  print *,'Valeurs min et max de X sur le maillage = ',minval(x),maxval(x)
+  print *,'Valeurs min et max de Z sur le maillage = ',minval(z),maxval(z)
   print *
 
+! afficher position de la source et des recepteurs
+  print *
+  print *,'Position (x,z) de la source'
+  print *
+  print *,'Source = ',xs,zs
+  print *
+  print *,'Position (x,z) des ',nrec,' recepteurs'
+  print *
+  do irec=1,nrec
+    print *,'Receiver ',irec,' = ',xrec(irec),zrec(irec)
+  enddo
+
 ! ***
-! *** generer un fichier 'GNUPLOT' pour le controle de la grille ***
+! *** generer un fichier Gnuplot pour le controle de la grille ***
 ! ***
 
   print *
-  print *,' Ecriture de la grille format GNUPLOT...'
+  print *,'Ecriture de la grille format Gnuplot...'
 
-  file1='gridfile.gnu'
+  open(unit=20,file='gridfile.gnu',status='unknown')
 
- open(unit=20,file=file1,status='unknown')
-
-! dessin de la topo de surface (splines)
-  do i=0,nx-1
-    write(20,15) sngl(absx(i)),sngl(top0(i))
-    write(20,15) sngl(absx(i+1)),sngl(top0(i+1))
-    write(20,10)
-  enddo
-
-! dessin de l'interface du milieu
-  if (interffile /= 'none') then
-  do i=0,nx-1
-    write(20,15) sngl(absx(i)),sngl(bot0(i))
-    write(20,15) sngl(absx(i+1)),sngl(bot0(i+1))
-    write(20,10)
-  enddo
-  endif
-
 ! dessin des lignes horizontales de la grille
-  print *, 'Ecriture lignes horizontales'
+  print *,'Ecriture lignes horizontales'
   istepx = 1
   if(ngnod == 4) then
     istepz = 1
@@ -597,7 +605,7 @@
   enddo
 
 ! dessin des lignes verticales de la grille
-  print *, 'Ecriture lignes verticales'
+  print *,'Ecriture lignes verticales'
   if(ngnod == 4) then
     istepx = 1
   else
@@ -615,13 +623,14 @@
   close(20)
 
 ! cree le script de dessin pour gnuplot
-  open(unit=20,file='plotgrid.gnu',status='unknown')
-  write(20,*) 'set term postscript landscape monochrome solid "Helvetica" 22'
-  write(20,*) 'set output "grille.ps"'
-  write(20,*) 'plot "gridfile.gnu" title "Macroblocs mesh" w l'
+  open(unit=20,file='plotgnu',status='unknown')
+  write(20,*) '#set term postscript landscape monochrome solid "Helvetica" 22'
+  write(20,*) '#set output "grille.ps"'
+  write(20,*) 'plot "gridfile.gnu" title "Macrobloc mesh" w l'
+  write(20,*) 'pause -1 "appuyez sur une touche"'
   close(20)
 
-  print *,' Fin ecriture de la grille format GNUPLOT'
+  print *,'Fin ecriture de la grille format Gnuplot'
   print *
 
 ! *** generation de la base de donnees
@@ -629,8 +638,8 @@
   open(unit=15,file='../SPECFEM90/DataBase',status='unknown')
 
   write(15,*) '#'
-  write(15,*) '# DataBase for SPECFEM2D version 5.0'
-  write(15,*) '# Dimitri Komatitsch, (c) University of Pau, France, May 2004'
+  write(15,*) '# DataBase for SPECFEM2D'
+  write(15,*) '# Dimitri Komatitsch, (c) University of Pau, France'
   write(15,*) '#'
 
   write(15,*) 'Titre simulation'
@@ -649,7 +658,7 @@
   write(15,*) ignuplot,interpol
 
   write(15,*) 'itaff icolor inumber'
-  write(15,*) itaff,0,0
+  write(15,*) itaff,1,0
 
   write(15,*) 'imeshvect imodelvect iboundvect cutvect isubsamp'
   write(15,*) imeshvect,imodelvect,iboundvect,cutvect,isubsamp
@@ -670,7 +679,7 @@
   write(15,*) nt,dt
 
   write(15,*) 'source'
-  write(15,*) '6 ',isource_type,xs-xmin,zs,f0,t0,factor,angle,' 0'
+  write(15,*) isource_type,xs-xmin,zs,f0,t0,factor,angle,' 0'
 
   write(15,*) 'Receiver positions:'
   do irec=1,nrec
@@ -703,12 +712,10 @@
   if(abshaut .and. absgauche) nelemabs = nelemabs - 1
   if(abshaut .and. absdroite) nelemabs = nelemabs - 1
 
-  netyp = 2
+  write(15,*) 'numat ngnod nspec iptsdisp ielemabs'
+  write(15,*) nbmodeles,ngnod,nspec,iptsdisp,nelemabs
 
-  write(15,*) 'netyp numat ngnod nspec iptsdisp ielemabs'
-  write(15,*) netyp,nbmodeles,ngnod,nspec,iptsdisp,nelemabs
-
-  write(15,*) 'Material sets (num 0 rho vp vs 0 0) or (num 2 rho c11 c13 c33 c44)'
+  write(15,*) 'Material sets (num 0 rho vp vs 0 0) or (num 1 rho c11 c13 c33 c44)'
   do i=1,nbmodeles
     write(15,*) i,icodemat(i),rho(i),cp(i),cs(i),aniso3(i),aniso4(i)
   enddo
@@ -759,17 +766,17 @@
   if(nelemabs > 0) then
   write(15,*) 'Liste des elements absorbants (haut bas gauche droite) :'
   inumabs = 0
-  do iz=1,nzread
-  do ix=1,nxread
+  do iz = 1,nzread
+  do ix = 1,nxread
     codehaut = .false.
     codebas = .false.
     codegauche = .false.
     codedroite = .false.
     inumelem = (iz-1)*nxread + ix
-    if(abshaut   .and. iz==nzread) codehaut = .true.
-    if(absbas    .and. iz== 1) codebas = .true.
-    if(absgauche .and. ix== 1) codegauche = .true.
-    if(absdroite .and. ix==nxread) codedroite = .true.
+    if(abshaut   .and. iz == nzread) codehaut = .true.
+    if(absbas    .and. iz == 1) codebas = .true.
+    if(absgauche .and. ix == 1) codegauche = .true.
+    if(absdroite .and. ix == nxread) codedroite = .true.
     if(codehaut .or. codebas .or. codegauche .or. codedroite) then
       inumabs = inumabs + 1
       write(15,*) inumabs,inumelem,codehaut,codebas,codegauche,codedroite
@@ -786,79 +793,54 @@
 
   end program meshfem2D
 
-! *****************
-! routines maillage
-! *****************
+! ********************
+! routines de maillage
+! ********************
 
 ! --- numero global du noeud
 
   integer function num(i,j,nx)
+
   implicit none
+
   integer i,j,nx
 
     num = j*(nx+1) + i + 1
 
   end function num
 
-! ------- definition des fonctions representant les interfaces -------
+! --- representation des interfaces par un spline
 
-!
-! --- bas du modele
-!
+  double precision function value_spline(x,xinterface,zinterface,coefs_interface,npoints_interface)
 
-  double precision function bottom(x)
   implicit none
-  double precision x
 
-    bottom = 0.d0
-
-  end function bottom
-
-!
-! --- representation interfaces par un spline
-!
-
-!--- spline
-
-  double precision function spl(x,xtopo,ztopo,coefs,ntopo)
-
-  implicit none
-  integer ntopo
+  integer npoints_interface
   double precision x,xp
-  double precision xtopo(ntopo),ztopo(ntopo)
-  double precision coefs(ntopo)
+  double precision, dimension(npoints_interface) :: xinterface,zinterface,coefs_interface
 
-  spl = 0.
+  value_spline = 0.d0
+
   xp = x
-  if (xp < xtopo(1)) xp = xtopo(1)
-  if (xp > xtopo(ntopo)) xp = xtopo(ntopo)
-  call splint(xtopo,ztopo,coefs,ntopo,xp,spl)
 
-  end function spl
+! si on sort du modele, prolonger par continuite
+  if(xp < xinterface(1)) xp = xinterface(1)
+  if(xp > xinterface(npoints_interface)) xp = xinterface(npoints_interface)
 
-! --- fonction de densification du maillage horizontal
+  call splint(xinterface,zinterface,coefs_interface,npoints_interface,xp,value_spline)
 
-  double precision function dens(ix,psi,xmin,xmax,nx)
+  end function value_spline
 
-  implicit none
-  integer ix,nx
-  double precision psi(0:nx)
-  double precision xmin,xmax
-
-  dens = xmin + dble(xmax-xmin)*psi(ix)
-
-  end function dens
-
 ! --------------------------------------
 
-! routine de calcul des coefs du spline (Numerical Recipes)
+! routine de calcul des coefs du spline (adapted from Numerical Recipes)
 
   subroutine spline(x,y,n,yp1,ypn,y2)
 
   implicit none
 
   integer n
-  double precision x(n),y(n),y2(n)
+  double precision, dimension(n) :: x,y,y2
   double precision, dimension(:), allocatable :: u
   double precision yp1,ypn
 
@@ -867,18 +849,21 @@
 
   allocate(u(n))
 
-  y2(1)=-0.5
-  u(1)=(3./(x(2)-x(1)))*((y(2)-y(1))/(x(2)-x(1))-yp1)
+  y2(1)=-0.5d0
+  u(1)=(3.d0/(x(2)-x(1)))*((y(2)-y(1))/(x(2)-x(1))-yp1)
+
   do i=2,n-1
     sig=(x(i)-x(i-1))/(x(i+1)-x(i-1))
-    p=sig*y2(i-1)+2.
-    y2(i)=(sig-1.)/p
-    u(i)=(6.*((y(i+1)-y(i))/(x(i+1)-x(i))-(y(i)-y(i-1)) &
+    p=sig*y2(i-1)+2.d0
+    y2(i)=(sig-1.d0)/p
+    u(i)=(6.d0*((y(i+1)-y(i))/(x(i+1)-x(i))-(y(i)-y(i-1)) &
                /(x(i)-x(i-1)))/(x(i+1)-x(i-1))-sig*u(i-1))/p
   enddo
-  qn=0.5
-  un=(3./(x(n)-x(n-1)))*(ypn-(y(n)-y(n-1))/(x(n)-x(n-1)))
-  y2(n)=(un-qn*u(n-1))/(qn*y2(n-1)+1.)
+
+  qn=0.5d0
+  un=(3.d0/(x(n)-x(n-1)))*(ypn-(y(n)-y(n-1))/(x(n)-x(n-1)))
+  y2(n)=(un-qn*u(n-1))/(qn*y2(n-1)+1.d0)
+
   do k=n-1,1,-1
     y2(k)=y2(k)*y2(k+1)+u(k)
   enddo
@@ -889,35 +874,38 @@
 
 ! --------------
 
-! routine d'evaluation du spline (Numerical Recipes)
+! routine d'evaluation du spline (adapted from Numerical Recipes)
 
-  SUBROUTINE SPLINT(XA,YA,Y2A,N,X,Y)
+  subroutine splint(XA,YA,Y2A,N,X,Y)
 
   implicit none
 
   integer n
-  double precision XA(N),YA(N),Y2A(N)
+  double precision, dimension(n) :: XA,YA,Y2A
   double precision x,y
 
   integer k,klo,khi
   double precision h,a,b
 
-  KLO=1
-  KHI=N
+  KLO = 1
+  KHI = N
+
   do while (KHI-KLO > 1)
     K=(KHI+KLO)/2
-    IF(XA(K) > X)THEN
+    if(XA(K) > X) then
       KHI=K
-    ELSE
+    else
       KLO=K
-    ENDIF
+    endif
   enddo
-  H=XA(KHI)-XA(KLO)
-  IF (H == 0.d0) stop 'Bad input in spline evaluation'
-  A=(XA(KHI)-X)/H
-  B=(X-XA(KLO))/H
 
-  Y=A*YA(KLO)+B*YA(KHI)+((A**3-A)*Y2A(KLO)+ (B**3-B)*Y2A(KHI))*(H**2)/6.d0
+  H = XA(KHI)-XA(KLO)
+  IF (H == 0.d0) stop 'bad input in spline evaluation'
 
-  end subroutine SPLINT
+  A = (XA(KHI)-X) / H
+  B = (X-XA(KLO)) / H
 
+  Y = A*YA(KLO) + B*YA(KHI) + ((A**3-A)*Y2A(KLO) + (B**3-B)*Y2A(KHI))*(H**2)/6.d0
+
+  end subroutine splint
+

Added: seismo/2D/SPECFEM2D/trunk/MAILLE90/old_topo_files/interf_paco.dat
===================================================================
--- seismo/2D/SPECFEM2D/trunk/MAILLE90/old_topo_files/interf_paco.dat	2004-12-14 15:22:01 UTC (rev 8432)
+++ seismo/2D/SPECFEM2D/trunk/MAILLE90/old_topo_files/interf_paco.dat	2007-12-07 23:45:56 UTC (rev 8433)
@@ -0,0 +1,212 @@
+   211
+   0.    4000.0000000000
+     33.333333333333    4000.0000000000
+     66.666666666667    4000.0000000000
+    100.000000000000    4000.0000000000
+     133.33333333333    4000.0000000000
+     166.66666666667    4000.0000000000
+     200.00000000000    4000.0000000000
+     233.33333333333    4000.0000000000
+     266.66666666667    4000.0000000000
+     300.00000000000    4000.0000000000
+     333.33333333333    4000.0000000000
+     366.66666666667    4000.0000000000
+     400.00000000000    4000.0000000000
+     433.33333333333    4000.0000000000
+     466.66666666667    4000.0000000000
+     500.00000000000    4000.0000000000
+     533.33333333333    4000.0000000000
+     566.66666666667    4000.0000000000
+     600.00000000000    4000.0000000000
+     633.33333333333    4000.0000000000
+     666.66666666667    4000.0000000000
+     700.00000000000    4000.0000000000
+     733.33333333333    4000.0000000000
+     766.66666666667    4000.0000000000
+     800.00000000000    4000.0000000000
+     833.33333333333    4000.0000000000
+     866.66666666667    4000.0000000000
+     900.00000000000    4000.0000000000
+     933.33333333333    4000.0000000000
+     966.66666666667    4000.0000000000
+    1000.00000000000    4000.0000000000
+     1033.3333333333    4000.0000000000
+     1066.6666666667    4000.0000000000
+     1100.0000000000    4000.0000000000
+     1133.3333333333    4000.0000000000
+     1166.6666666667    4000.0000000000
+     1200.0000000000    4000.0000000000
+     1233.3333333333    4000.0000000000
+     1266.6666666667    4000.0000000000
+     1300.0000000000    4000.0000000000
+     1333.3333333333    4000.0000000000
+     1366.6666666667    4000.0000000000
+     1400.0000000000    4000.0000000000
+     1433.3333333333    4000.0000000000
+     1466.6666666667    4000.0000000000
+     1500.0000000000    4000.0000000000
+     1533.3333333333    4000.0000000000
+     1566.6666666667    4000.0000000000
+     1600.0000000000    4000.0000000000
+     1633.3333333333    4000.0000000000
+     1666.6666666667    4000.0000000000
+     1700.0000000000    4000.0000000000
+     1733.3333333333    4000.0000000000
+     1766.6666666667    4000.0000000000
+     1800.0000000000    4000.0000000000
+     1833.3333333333    4000.0000000000
+     1866.6666666667    4000.0000000000
+     1900.0000000000    4000.0000000000
+     1933.3333333333    4000.0000000000
+     1966.6666666667    4000.0000000000
+     2000.0000000000    4000.0000000000
+     2033.3333333333    4000.0000000000
+     2066.6666666667    4000.0000000000
+     2100.0000000000    4000.0000000000
+     2133.3333333333    4000.0000000000
+     2166.6666666667    4000.0000000000
+     2200.0000000000    4000.0000000000
+     2233.3333333333    4000.0000000000
+     2266.6666666667    4000.0000000000
+     2300.0000000000    4000.0000000000
+     2333.3333333333    4000.0000000000
+     2366.6666666667    4000.0000000000
+     2400.0000000000    4000.0000000000
+     2433.3333333333    4000.0000000000
+     2466.6666666667    4000.0000000000
+     2500.0000000000    4000.0000000000
+     2533.3333333333    3998.6304738452
+     2566.6666666667    3994.5369001959
+     2600.0000000000    3987.7641291015
+     2633.3333333333    3978.3863644593
+     2666.6666666667    3966.5063510209
+     2700.0000000000    3952.2542486992
+     2733.3333333333    3935.7862065095
+     2766.6666666667    3917.2826517676
+     2800.0000000000    3896.9463132909
+     2833.3333333333    3875.0000002591
+     2866.6666666667    3851.6841610696
+     2900.0000000000    3827.2542489351
+     2933.3333333333    3801.9779230848
+     2966.6666666667    3776.1321162334
+     3000.0000000000    3750.0000004487
+     3033.3333333333    3723.8678846591
+     3066.6666666667    3698.0220777930
+     3100.0000000000    3672.7457519184
+     3133.3333333333    3648.3158397503
+     3166.6666666667    3625.0000005181
+     3200.0000000000    3603.0536874351
+     3233.3333333333    3582.7173488994
+     3266.6666666667    3564.2137940910
+     3300.0000000000    3547.7457518283
+     3333.3333333333    3533.4936494278
+     3366.6666666667    3521.6136359057
+     3400.0000000000    3512.2358711758
+     3433.3333333333    3505.4630999907
+     3466.6666666667    3501.3695262486
+     3500.0000000000    3500.0000000000
+     3533.3333333333    3501.3695260610
+     3566.6666666667    3505.4630996175
+     3600.0000000000    3512.2358706212
+     3633.3333333333    3521.6136351757
+     3666.6666666667    3533.4936485304
+     3700.0000000000    3547.7457507733
+     3733.3333333333    3564.2137928900
+     3766.6666666667    3582.7173475655
+     3800.0000000000    3603.0536859830
+     3833.3333333333    3624.9999989637
+     3866.6666666667    3648.3158381106
+     3900.0000000000    3672.7457502113
+     3933.3333333333    3698.0220760373
+     3966.6666666667    3723.8678828740
+     4000.0000000000    3749.9999986538
+     4033.3333333333    3776.1321144484
+     4066.6666666667    3801.9779213292
+     4100.0000000000    3827.2542472281
+     4133.3333333333    3851.6841594298
+     4166.6666666667    3874.9999987046
+     4200.0000000000    3896.9463118388
+     4233.3333333333    3917.2826504337
+     4266.6666666667    3935.7862053084
+     4300.0000000000    3952.2542476442
+     4333.3333333333    3966.5063501234
+     4366.6666666667    3978.3863637293
+     4400.0000000000    3987.7641285469
+     4433.3333333333    3994.5368998227
+     4466.6666666667    3998.6304736576
+     4500.0000000000    4000.0000000000
+     4533.3333333333    4000.0000000000
+     4566.6666666667    4000.0000000000
+     4600.0000000000    4000.0000000000
+     4633.3333333333    4000.0000000000
+     4666.6666666667    4000.0000000000
+     4700.0000000000    4000.0000000000
+     4733.3333333333    4000.0000000000
+     4766.6666666667    4000.0000000000
+     4800.0000000000    4000.0000000000
+     4833.3333333333    4000.0000000000
+     4866.6666666667    4000.0000000000
+     4900.0000000000    4000.0000000000
+     4933.3333333333    4000.0000000000
+     4966.6666666667    4000.0000000000
+     5000.0000000000    4000.0000000000
+     5033.3333333333    4000.0000000000
+     5066.6666666667    4000.0000000000
+     5100.0000000000    4000.0000000000
+     5133.3333333333    4000.0000000000
+     5166.6666666667    4000.0000000000
+     5200.0000000000    4000.0000000000
+     5233.3333333333    4000.0000000000
+     5266.6666666667    4000.0000000000
+     5300.0000000000    4000.0000000000
+     5333.3333333333    4000.0000000000
+     5366.6666666667    4000.0000000000
+     5400.0000000000    4000.0000000000
+     5433.3333333333    4000.0000000000
+     5466.6666666667    4000.0000000000
+     5500.0000000000    4000.0000000000
+     5533.3333333333    4000.0000000000
+     5566.6666666667    4000.0000000000
+     5600.0000000000    4000.0000000000
+     5633.3333333333    4000.0000000000
+     5666.6666666667    4000.0000000000
+     5700.0000000000    4000.0000000000
+     5733.3333333333    4000.0000000000
+     5766.6666666667    4000.0000000000
+     5800.0000000000    4000.0000000000
+     5833.3333333333    4000.0000000000
+     5866.6666666667    4000.0000000000
+     5900.0000000000    4000.0000000000
+     5933.3333333333    4000.0000000000
+     5966.6666666667    4000.0000000000
+     6000.0000000000    4000.0000000000
+     6033.3333333333    4000.0000000000
+     6066.6666666667    4000.0000000000
+     6100.0000000000    4000.0000000000
+     6133.3333333333    4000.0000000000
+     6166.6666666667    4000.0000000000
+     6200.0000000000    4000.0000000000
+     6233.3333333333    4000.0000000000
+     6266.6666666667    4000.0000000000
+     6300.0000000000    4000.0000000000
+     6333.3333333333    4000.0000000000
+     6366.6666666667    4000.0000000000
+     6400.0000000000    4000.0000000000
+     6433.3333333333    4000.0000000000
+     6466.6666666667    4000.0000000000
+     6500.0000000000    4000.0000000000
+     6533.3333333333    4000.0000000000
+     6566.6666666667    4000.0000000000
+     6600.0000000000    4000.0000000000
+     6633.3333333333    4000.0000000000
+     6666.6666666667    4000.0000000000
+     6700.0000000000    4000.0000000000
+     6733.3333333333    4000.0000000000
+     6766.6666666667    4000.0000000000
+     6800.0000000000    4000.0000000000
+     6833.3333333333    4000.0000000000
+     6866.6666666667    4000.0000000000
+     6900.0000000000    4000.0000000000
+     6933.3333333333    4000.0000000000
+     6966.6666666667    4000.0000000000
+     7000.0000000000    4000.0000000000

Added: seismo/2D/SPECFEM2D/trunk/MAILLE90/old_topo_files/topo_cours_M2_UPPA_curved.dat
===================================================================
--- seismo/2D/SPECFEM2D/trunk/MAILLE90/old_topo_files/topo_cours_M2_UPPA_curved.dat	2004-12-14 15:22:01 UTC (rev 8432)
+++ seismo/2D/SPECFEM2D/trunk/MAILLE90/old_topo_files/topo_cours_M2_UPPA_curved.dat	2007-12-07 23:45:56 UTC (rev 8433)
@@ -0,0 +1,10 @@
+9
+0 3000
+500 3000
+1000 3100
+1500 3400
+2000 3300
+2500 3200
+3000 3100
+3500 3000
+5000 3000

Added: seismo/2D/SPECFEM2D/trunk/MAILLE90/old_topo_files/topo_rouen_aniso.dat
===================================================================
--- seismo/2D/SPECFEM2D/trunk/MAILLE90/old_topo_files/topo_rouen_aniso.dat	2004-12-14 15:22:01 UTC (rev 8432)
+++ seismo/2D/SPECFEM2D/trunk/MAILLE90/old_topo_files/topo_rouen_aniso.dat	2007-12-07 23:45:56 UTC (rev 8433)
@@ -0,0 +1,3 @@
+2
+0 0.06d0
+1 0.06d0

Added: seismo/2D/SPECFEM2D/trunk/MAILLE90/old_topo_files/topoarticle.dat
===================================================================
--- seismo/2D/SPECFEM2D/trunk/MAILLE90/old_topo_files/topoarticle.dat	2004-12-14 15:22:01 UTC (rev 8432)
+++ seismo/2D/SPECFEM2D/trunk/MAILLE90/old_topo_files/topoarticle.dat	2007-12-07 23:45:56 UTC (rev 8433)
@@ -0,0 +1,210 @@
+ 209
+   0.    1741.09
+     12.0192    1743.55
+     24.0385    1746.09
+     36.0577    1748.72
+     48.0769    1751.41
+     60.0962    1754.16
+     72.1154    1756.95
+     84.1346    1759.79
+     96.1538    1762.66
+     108.173    1765.56
+     120.192    1768.47
+     132.212    1771.38
+     144.231    1774.30
+     156.250    1777.21
+     168.269    1780.11
+     180.288    1782.99
+     192.308    1785.85
+     204.327    1788.68
+     216.346    1791.46
+     228.365    1794.20
+     240.385    1796.89
+     252.404    1799.52
+     264.423    1802.09
+     276.442    1804.60
+     288.462    1807.02
+     300.481    1809.37
+     312.500    1811.63
+     324.519    1813.80
+     336.538    1815.90
+     348.558    1817.91
+     360.577    1819.86
+     372.596    1821.73
+     384.615    1823.53
+     396.635    1825.27
+     408.654    1826.95
+     420.673    1828.58
+     432.692    1830.19
+     444.712    1831.78
+     456.731    1833.38
+     468.750    1835.00
+     480.769    1836.67
+     492.788    1838.40
+     504.808    1840.21
+     516.827    1842.10
+     528.846    1844.02
+     540.865    1845.95
+     552.885    1847.86
+     564.904    1849.70
+     576.923    1851.44
+     588.942    1853.05
+     600.962    1854.49
+     612.981    1855.74
+     625.000    1856.81
+     637.019    1857.72
+     649.038    1858.49
+     661.058    1859.13
+     673.077    1859.68
+     685.096    1860.15
+     697.115    1860.56
+     709.135    1860.92
+     721.154    1861.22
+     733.173    1861.42
+     745.192    1861.49
+     757.212    1861.39
+     769.231    1861.09
+     781.250    1860.56
+     793.269    1859.75
+     805.288    1858.64
+     817.308    1857.17
+     829.327    1855.28
+     841.346    1852.89
+     853.365    1849.93
+     865.385    1846.34
+     877.404    1842.04
+     889.423    1836.97
+     901.442    1831.04
+     913.462    1824.28
+     925.481    1816.85
+     937.500    1809.00
+     949.519    1800.95
+     961.538    1792.92
+     973.558    1785.14
+     985.577    1777.84
+     997.596    1771.25
+     1009.62    1765.55
+     1021.63    1760.76
+     1033.65    1756.75
+     1045.67    1753.44
+     1057.69    1750.72
+     1069.71    1748.49
+     1081.73    1746.65
+     1093.75    1745.09
+     1105.77    1743.72
+     1117.79    1742.49
+     1129.81    1741.37
+     1141.83    1740.35
+     1153.85    1739.42
+     1165.87    1738.57
+     1177.88    1737.79
+     1189.90    1737.06
+     1201.92    1736.37
+     1213.94    1735.72
+     1225.96    1735.14
+     1237.98    1734.66
+     1250.00    1734.30
+     1262.02    1734.11
+     1274.04    1734.11
+     1286.06    1734.33
+     1298.08    1734.80
+     1310.10    1735.55
+     1322.12    1736.55
+     1334.13    1737.74
+     1346.15    1739.08
+     1358.17    1740.50
+     1370.19    1741.96
+     1382.21    1743.40
+     1394.23    1744.77
+     1406.25    1746.02
+     1418.27    1747.13
+     1430.29    1748.10
+     1442.31    1748.94
+     1454.33    1749.64
+     1466.35    1750.22
+     1478.37    1750.68
+     1490.38    1751.02
+     1502.40    1751.25
+     1514.42    1751.41
+     1526.44    1751.68
+     1538.46    1752.23
+     1550.48    1753.22
+     1562.50    1754.84
+     1574.52    1757.27
+     1586.54    1760.67
+     1598.56    1765.22
+     1610.58    1771.04
+     1622.60    1777.93
+     1634.62    1785.56
+     1646.63    1793.60
+     1658.65    1801.74
+     1670.67    1809.64
+     1682.69    1817.00
+     1694.71    1823.48
+     1706.73    1828.77
+     1718.75    1832.75
+     1730.77    1835.44
+     1742.79    1836.88
+     1754.81    1837.11
+     1766.83    1836.17
+     1778.85    1834.08
+     1790.87    1830.90
+     1802.88    1826.64
+     1814.90    1821.41
+     1826.92    1815.42
+     1838.94    1808.87
+     1850.96    1802.01
+     1862.98    1795.03
+     1875.00    1788.17
+     1887.02    1781.64
+     1899.04    1775.67
+     1911.06    1770.45
+     1923.08    1765.99
+     1935.10    1762.29
+     1947.12    1759.31
+     1959.13    1757.05
+     1971.15    1755.47
+     1983.17    1754.55
+     1995.19    1754.28
+     2007.21    1754.63
+     2019.23    1755.57
+     2031.25    1757.03
+     2043.27    1758.96
+     2055.29    1761.33
+     2067.31    1764.06
+     2079.33    1767.11
+     2091.35    1770.43
+     2103.37    1773.97
+     2115.38    1777.65
+     2127.40    1781.39
+     2139.42    1785.08
+     2151.44    1788.64
+     2163.46    1791.97
+     2175.48    1794.97
+     2187.50    1797.54
+     2199.52    1799.60
+     2211.54    1801.07
+     2223.56    1802.01
+     2235.58    1802.51
+     2247.60    1802.66
+     2259.62    1802.56
+     2271.63    1802.29
+     2283.65    1801.94
+     2295.67    1801.62
+     2307.69    1801.39
+     2319.71    1801.30
+     2331.73    1801.30
+     2343.75    1801.39
+     2355.77    1801.54
+     2367.79    1801.73
+     2379.81    1801.93
+     2391.83    1802.13
+     2403.85    1802.30
+     2415.87    1802.43
+     2427.88    1802.51
+     2439.90    1802.55
+     2451.92    1802.56
+     2463.94    1802.52
+     2475.96    1802.45
+     2487.98    1802.35
+     2500.00    1802.20

Deleted: seismo/2D/SPECFEM2D/trunk/MAILLE90/profilx.dat
===================================================================
--- seismo/2D/SPECFEM2D/trunk/MAILLE90/profilx.dat	2004-12-14 15:22:01 UTC (rev 8432)
+++ seismo/2D/SPECFEM2D/trunk/MAILLE90/profilx.dat	2007-12-07 23:45:56 UTC (rev 8433)
@@ -1,215 +0,0 @@
--4.030060
--4.031510
--4.027550
--4.022690
--4.017840
--4.011910
--3.999230
--3.972270
--3.940420
--3.909460
--3.885800
--3.868120
--3.853440
--3.840630
--3.830040
--3.818850
--3.804520
--3.790030
--3.776320
--3.765340
--3.755920
--3.746830
--3.735490
--3.721820
--3.703930
--3.679520
--3.650700
--3.620200
--3.591410
--3.568040
--3.547570
--3.528300
--3.507120
--3.485880
--3.468120
--3.446460
--3.420720
--3.390030
--3.353750
--3.313270
--3.264390
--3.214030
--3.179940
--3.155080
--3.132110
--3.106450
--3.075090
--3.038830
--2.991240
--2.916530
--2.815720
--2.748110
--2.712470
--2.670830
--2.618150
--2.555080
--2.484160
--2.405320
--2.321900
--2.229520
--2.133560
--2.070170
--2.027820
--1.974430
--1.907660
--1.851900
--1.790190
--1.724810
--1.662320
--1.613210
--1.567500
--1.532530
--1.499660
--1.420860
--1.326870
--1.238010
--1.149630
--1.055670
--0.951941
--0.830088
--0.695442
--0.535123
--0.303787
--0.066964
-0.137585
-0.292947
-0.427345
-0.504895
-0.548948
-0.611893
-0.669241
-0.731076
-0.784069
-0.862715
-0.908088
-0.989800
-1.074750
-1.125830
-1.155250
-1.142880
-1.068660
-0.966157
-1.106490
-1.229600
-1.452710
-1.620460
-1.795150
-1.775980
-1.694660
-1.678880
-1.661540
-1.600500
-1.536210
-1.454460
-1.203100
-0.999440
-0.972028
-0.935924
-0.995477
-1.063760
-1.098300
-1.056350
-0.986176
-0.855505
-0.689433
-0.529934
-0.397784
-0.317342
-0.222705
-0.117320
-0.028807
--0.040407
--0.118587
--0.326905
--0.574415
--0.741248
--0.887590
--1.010380
--1.087830
--1.108840
--1.160430
--1.209690
--1.264010
--1.319100
--1.348980
--1.404960
--1.485740
--1.575370
--1.662450
--1.738790
--1.798140
--1.850290
--1.905520
--1.951300
--1.995390
--2.063230
--2.151310
--2.217030
--2.278450
--2.333810
--2.384270
--2.421330
--2.443130
--2.465170
--2.567050
--2.690890
--2.770720
--2.838250
--2.898530
--2.956700
--3.026270
--3.102020
--3.173380
--3.233160
--3.285530
--3.338610
--3.390250
--3.433800
--3.482320
--3.537380
--3.588740
--3.628460
--3.662240
--3.688980
--3.709740
--3.732810
--3.750350
--3.772050
--3.797530
--3.827180
--3.859440
--3.890500
--3.918430
--3.944920
--3.968160
--3.987850
--4.002890
--4.012180
--4.014820
--4.016240
--4.021080
--4.025060
--4.027960
--4.030390
--4.031410
--4.032690
--4.034570
--4.035590
--4.036420
--4.037270
--4.037550
--4.016490
--4.008760
--4.007300
--4.006200

Deleted: seismo/2D/SPECFEM2D/trunk/MAILLE90/profily.dat
===================================================================
--- seismo/2D/SPECFEM2D/trunk/MAILLE90/profily.dat	2004-12-14 15:22:01 UTC (rev 8432)
+++ seismo/2D/SPECFEM2D/trunk/MAILLE90/profily.dat	2007-12-07 23:45:56 UTC (rev 8433)
@@ -1,209 +0,0 @@
--4.042630
--4.033770
--4.024570
--4.031010
--4.027090
--4.011370
--3.989330
--3.963060
--3.929730
--3.893400
--3.854910
--3.807280
--3.764700
--3.718700
--3.667760
--3.619570
--3.579390
--3.529760
--3.477630
--3.427720
--3.383820
--3.346060
--3.312440
--3.281950
--3.252910
--3.222660
--3.189850
--3.153580
--3.113840
--3.073560
--3.036520
--2.994700
--2.946600
--2.906230
--2.856750
--2.801230
--2.751650
--2.709890
--2.593630
--2.501410
--2.339970
--2.219030
--2.157610
--2.126910
--2.081700
--2.033360
--1.993200
--1.959890
--1.920950
--1.885440
--1.849060
--1.795820
--1.723950
--1.624180
--1.500000
--1.379030
--1.266730
--1.161710
--1.064370
--0.991573
--0.971158
--0.935858
--0.900854
--0.874644
--0.842509
--0.808727
--0.789044
--0.797708
--0.768248
--0.756078
--0.738104
--0.756372
--0.775592
--0.795957
--0.787438
--0.705082
--0.595518
--0.438898
--0.292028
--0.145967
--0.062730
-0.042711
-0.270900
-0.520289
-0.726283
-0.910875
-1.135880
-1.376630
-1.592740
-1.757490
-1.894770
-2.043840
-2.106480
-2.091720
-1.818020
-1.200320
-0.943632
-0.864765
-0.849034
-0.953888
-1.012090
-1.158530
-1.759700
-1.777570
-1.252260
-1.043510
-1.229600
-1.496690
-1.515210
-1.522500
-1.522040
-1.505200
-1.497920
-1.479140
-1.447090
-1.430400
-1.462890
-1.560170
-1.666450
-1.695220
-1.675090
-1.617180
-1.532760
-1.426300
-1.375350
-1.198940
-0.876019
-0.921095
-0.998579
-0.847144
-0.712094
-0.580680
-0.459470
-0.335168
-0.179500
-0.065837
--0.008686
--0.125625
--0.407411
--0.779719
--1.215620
--1.613050
--1.707880
--1.708050
--1.734230
--1.753980
--1.773940
--1.796890
--1.822940
--1.854070
--1.880990
--1.901660
--1.965350
--1.995410
--1.941710
--1.898150
--1.983790
--2.167880
--2.287980
--2.395880
--2.460820
--2.487850
--2.504990
--2.499650
--2.511300
--2.511300
--2.551220
--2.608210
--2.644140
--2.692670
--2.720570
--2.735400
--2.747100
--2.788430
--2.834740
--2.886570
--2.951890
--3.016690
--3.075140
--3.116640
--3.156450
--3.207550
--3.258550
--3.308610
--3.356020
--3.398640
--3.437490
--3.477580
--3.513470
--3.544880
--3.568640
--3.600160
--3.632250
--3.665310
--3.697350
--3.729730
--3.763660
--3.795090
--3.823330
--3.854270
--3.883290
--3.908390
--3.931510
--3.955920
--3.977590
--3.992900
--4.004930
--4.014750
--4.025040

Deleted: seismo/2D/SPECFEM2D/trunk/MAILLE90/topo_cours_M2_UPPA_curved.dat
===================================================================
--- seismo/2D/SPECFEM2D/trunk/MAILLE90/topo_cours_M2_UPPA_curved.dat	2004-12-14 15:22:01 UTC (rev 8432)
+++ seismo/2D/SPECFEM2D/trunk/MAILLE90/topo_cours_M2_UPPA_curved.dat	2007-12-07 23:45:56 UTC (rev 8433)
@@ -1,10 +0,0 @@
-9
-0 3000
-500 3000
-1000 3100
-1500 3400
-2000 3300
-2500 3200
-3000 3100
-3500 3000
-5000 3000

Deleted: seismo/2D/SPECFEM2D/trunk/MAILLE90/topo_rouen_aniso.dat
===================================================================
--- seismo/2D/SPECFEM2D/trunk/MAILLE90/topo_rouen_aniso.dat	2004-12-14 15:22:01 UTC (rev 8432)
+++ seismo/2D/SPECFEM2D/trunk/MAILLE90/topo_rouen_aniso.dat	2007-12-07 23:45:56 UTC (rev 8433)
@@ -1,3 +0,0 @@
-2
-0 0.06d0
-1 0.06d0

Deleted: seismo/2D/SPECFEM2D/trunk/MAILLE90/topoarticle.dat
===================================================================
--- seismo/2D/SPECFEM2D/trunk/MAILLE90/topoarticle.dat	2004-12-14 15:22:01 UTC (rev 8432)
+++ seismo/2D/SPECFEM2D/trunk/MAILLE90/topoarticle.dat	2007-12-07 23:45:56 UTC (rev 8433)
@@ -1,210 +0,0 @@
- 209
-   0.    1741.09
-     12.0192    1743.55
-     24.0385    1746.09
-     36.0577    1748.72
-     48.0769    1751.41
-     60.0962    1754.16
-     72.1154    1756.95
-     84.1346    1759.79
-     96.1538    1762.66
-     108.173    1765.56
-     120.192    1768.47
-     132.212    1771.38
-     144.231    1774.30
-     156.250    1777.21
-     168.269    1780.11
-     180.288    1782.99
-     192.308    1785.85
-     204.327    1788.68
-     216.346    1791.46
-     228.365    1794.20
-     240.385    1796.89
-     252.404    1799.52
-     264.423    1802.09
-     276.442    1804.60
-     288.462    1807.02
-     300.481    1809.37
-     312.500    1811.63
-     324.519    1813.80
-     336.538    1815.90
-     348.558    1817.91
-     360.577    1819.86
-     372.596    1821.73
-     384.615    1823.53
-     396.635    1825.27
-     408.654    1826.95
-     420.673    1828.58
-     432.692    1830.19
-     444.712    1831.78
-     456.731    1833.38
-     468.750    1835.00
-     480.769    1836.67
-     492.788    1838.40
-     504.808    1840.21
-     516.827    1842.10
-     528.846    1844.02
-     540.865    1845.95
-     552.885    1847.86
-     564.904    1849.70
-     576.923    1851.44
-     588.942    1853.05
-     600.962    1854.49
-     612.981    1855.74
-     625.000    1856.81
-     637.019    1857.72
-     649.038    1858.49
-     661.058    1859.13
-     673.077    1859.68
-     685.096    1860.15
-     697.115    1860.56
-     709.135    1860.92
-     721.154    1861.22
-     733.173    1861.42
-     745.192    1861.49
-     757.212    1861.39
-     769.231    1861.09
-     781.250    1860.56
-     793.269    1859.75
-     805.288    1858.64
-     817.308    1857.17
-     829.327    1855.28
-     841.346    1852.89
-     853.365    1849.93
-     865.385    1846.34
-     877.404    1842.04
-     889.423    1836.97
-     901.442    1831.04
-     913.462    1824.28
-     925.481    1816.85
-     937.500    1809.00
-     949.519    1800.95
-     961.538    1792.92
-     973.558    1785.14
-     985.577    1777.84
-     997.596    1771.25
-     1009.62    1765.55
-     1021.63    1760.76
-     1033.65    1756.75
-     1045.67    1753.44
-     1057.69    1750.72
-     1069.71    1748.49
-     1081.73    1746.65
-     1093.75    1745.09
-     1105.77    1743.72
-     1117.79    1742.49
-     1129.81    1741.37
-     1141.83    1740.35
-     1153.85    1739.42
-     1165.87    1738.57
-     1177.88    1737.79
-     1189.90    1737.06
-     1201.92    1736.37
-     1213.94    1735.72
-     1225.96    1735.14
-     1237.98    1734.66
-     1250.00    1734.30
-     1262.02    1734.11
-     1274.04    1734.11
-     1286.06    1734.33
-     1298.08    1734.80
-     1310.10    1735.55
-     1322.12    1736.55
-     1334.13    1737.74
-     1346.15    1739.08
-     1358.17    1740.50
-     1370.19    1741.96
-     1382.21    1743.40
-     1394.23    1744.77
-     1406.25    1746.02
-     1418.27    1747.13
-     1430.29    1748.10
-     1442.31    1748.94
-     1454.33    1749.64
-     1466.35    1750.22
-     1478.37    1750.68
-     1490.38    1751.02
-     1502.40    1751.25
-     1514.42    1751.41
-     1526.44    1751.68
-     1538.46    1752.23
-     1550.48    1753.22
-     1562.50    1754.84
-     1574.52    1757.27
-     1586.54    1760.67
-     1598.56    1765.22
-     1610.58    1771.04
-     1622.60    1777.93
-     1634.62    1785.56
-     1646.63    1793.60
-     1658.65    1801.74
-     1670.67    1809.64
-     1682.69    1817.00
-     1694.71    1823.48
-     1706.73    1828.77
-     1718.75    1832.75
-     1730.77    1835.44
-     1742.79    1836.88
-     1754.81    1837.11
-     1766.83    1836.17
-     1778.85    1834.08
-     1790.87    1830.90
-     1802.88    1826.64
-     1814.90    1821.41
-     1826.92    1815.42
-     1838.94    1808.87
-     1850.96    1802.01
-     1862.98    1795.03
-     1875.00    1788.17
-     1887.02    1781.64
-     1899.04    1775.67
-     1911.06    1770.45
-     1923.08    1765.99
-     1935.10    1762.29
-     1947.12    1759.31
-     1959.13    1757.05
-     1971.15    1755.47
-     1983.17    1754.55
-     1995.19    1754.28
-     2007.21    1754.63
-     2019.23    1755.57
-     2031.25    1757.03
-     2043.27    1758.96
-     2055.29    1761.33
-     2067.31    1764.06
-     2079.33    1767.11
-     2091.35    1770.43
-     2103.37    1773.97
-     2115.38    1777.65
-     2127.40    1781.39
-     2139.42    1785.08
-     2151.44    1788.64
-     2163.46    1791.97
-     2175.48    1794.97
-     2187.50    1797.54
-     2199.52    1799.60
-     2211.54    1801.07
-     2223.56    1802.01
-     2235.58    1802.51
-     2247.60    1802.66
-     2259.62    1802.56
-     2271.63    1802.29
-     2283.65    1801.94
-     2295.67    1801.62
-     2307.69    1801.39
-     2319.71    1801.30
-     2331.73    1801.30
-     2343.75    1801.39
-     2355.77    1801.54
-     2367.79    1801.73
-     2379.81    1801.93
-     2391.83    1802.13
-     2403.85    1802.30
-     2415.87    1802.43
-     2427.88    1802.51
-     2439.90    1802.55
-     2451.92    1802.56
-     2463.94    1802.52
-     2475.96    1802.45
-     2487.98    1802.35
-     2500.00    1802.20

Modified: seismo/2D/SPECFEM2D/trunk/SPECFEM90/Makefile
===================================================================
--- seismo/2D/SPECFEM2D/trunk/SPECFEM90/Makefile	2004-12-14 15:22:01 UTC (rev 8432)
+++ seismo/2D/SPECFEM2D/trunk/SPECFEM90/Makefile	2007-12-07 23:45:56 UTC (rev 8433)
@@ -1,7 +1,7 @@
 #
-# Makefile for SPECFEM2D version 5.0
+# Makefile for SPECFEM2D version 5.1
 #
-# Dimitri Komatitsch, Universite de Pau et des Pays de l'Adour, May 2004
+# Dimitri Komatitsch, Universite de Pau et des Pays de l'Adour, December 2004
 # 
 SHELL=/bin/sh
 
@@ -13,7 +13,7 @@
 
 # Intel Linux
 F90 = ifort
-#FLAGS=-O0 -implicitnone -warn stderrors -warn truncated_source -warn argument_checking -warn unused -warn declarations -std95 -check bounds -C
+#FLAGS=-O0 -implicitnone -warn stderrors -warn truncated_source -warn argument_checking -warn unused -warn declarations -std95 -check bounds -CB
 FLAGS=-O2 -implicitnone -warn stderrors -warn truncated_source -warn argument_checking -warn unused -warn declarations -std95 -check nobounds
 
 #

Modified: seismo/2D/SPECFEM2D/trunk/SPECFEM90/checkgrid.f90
===================================================================
--- seismo/2D/SPECFEM2D/trunk/SPECFEM90/checkgrid.f90	2004-12-14 15:22:01 UTC (rev 8432)
+++ seismo/2D/SPECFEM2D/trunk/SPECFEM90/checkgrid.f90	2007-12-07 23:45:56 UTC (rev 8433)
@@ -1,13 +1,13 @@
 
 !========================================================================
 !
-!                   S P E C F E M 2 D  Version 5.0
+!                   S P E C F E M 2 D  Version 5.1
 !                   ------------------------------
 !
 !                         Dimitri Komatitsch
 !          Universite de Pau et des Pays de l'Adour, France
 !
-!                          (c) May 2004
+!                          (c) December 2004
 !
 !========================================================================
 

Modified: seismo/2D/SPECFEM2D/trunk/SPECFEM90/compute_gradient_attenuation.f90
===================================================================
--- seismo/2D/SPECFEM2D/trunk/SPECFEM90/compute_gradient_attenuation.f90	2004-12-14 15:22:01 UTC (rev 8432)
+++ seismo/2D/SPECFEM2D/trunk/SPECFEM90/compute_gradient_attenuation.f90	2007-12-07 23:45:56 UTC (rev 8433)
@@ -1,13 +1,13 @@
 
 !========================================================================
 !
-!                   S P E C F E M 2 D  Version 5.0
+!                   S P E C F E M 2 D  Version 5.1
 !                   ------------------------------
 !
 !                         Dimitri Komatitsch
 !          Universite de Pau et des Pays de l'Adour, France
 !
-!                          (c) May 2004
+!                          (c) December 2004
 !
 !========================================================================
 

Modified: seismo/2D/SPECFEM2D/trunk/SPECFEM90/createnum_fast.f90
===================================================================
--- seismo/2D/SPECFEM2D/trunk/SPECFEM90/createnum_fast.f90	2004-12-14 15:22:01 UTC (rev 8432)
+++ seismo/2D/SPECFEM2D/trunk/SPECFEM90/createnum_fast.f90	2007-12-07 23:45:56 UTC (rev 8433)
@@ -1,13 +1,13 @@
 
 !========================================================================
 !
-!                   S P E C F E M 2 D  Version 5.0
+!                   S P E C F E M 2 D  Version 5.1
 !                   ------------------------------
 !
 !                         Dimitri Komatitsch
 !          Universite de Pau et des Pays de l'Adour, France
 !
-!                          (c) May 2004
+!                          (c) December 2004
 !
 !========================================================================
 
@@ -43,8 +43,8 @@
   print *,'Generating global mesh numbering (fast version)...'
   print *
 
-  nxyz   = NGLLX*NGLLZ
-  ntot   = nxyz*nspec
+  nxyz = NGLLX*NGLLZ
+  ntot = nxyz*nspec
 
   allocate(loc(ntot))
   allocate(ind(ntot))
@@ -104,15 +104,15 @@
   ymaxval=-HUGEVAL
   ieoff=nxyz*(ie-1)
   do ilocnum=1,nxyz
-    xmaxval=dmax1(xp(ieoff+ilocnum),xmaxval)
-    xminval=dmin1(xp(ieoff+ilocnum),xminval)
-    ymaxval=dmax1(yp(ieoff+ilocnum),ymaxval)
-    yminval=dmin1(yp(ieoff+ilocnum),yminval)
+    xmaxval=max(xp(ieoff+ilocnum),xmaxval)
+    xminval=min(xp(ieoff+ilocnum),xminval)
+    ymaxval=max(yp(ieoff+ilocnum),ymaxval)
+    yminval=min(yp(ieoff+ilocnum),yminval)
   enddo
 
 ! compute the minimum typical "size" of an element in the mesh
-  xtypdist = dmin1(xtypdist,xmaxval-xminval)
-  xtypdist = dmin1(xtypdist,ymaxval-yminval)
+  xtypdist = min(xtypdist,xmaxval-xminval)
+  xtypdist = min(xtypdist,ymaxval-yminval)
 
   enddo
 
@@ -128,10 +128,10 @@
 !  Sort within each segment
    ioff=1
    do iseg=1,nseg
-      if  (j == 1) then
-         call rank (xp(ioff),ind,ninseg(iseg))
+      if(j == 1) then
+        call rank (xp(ioff),ind,ninseg(iseg))
       else
-         call rank (yp(ioff),ind,ninseg(iseg))
+        call rank (yp(ioff),ind,ninseg(iseg))
       endif
       call swap(xp(ioff),work,ind,ninseg(iseg))
       call swap(yp(ioff),work,ind,ninseg(iseg))
@@ -211,6 +211,7 @@
 ! Use Heap Sort (p 233 Numerical Recipes)
 !
   implicit none
+
   integer N
   double precision A(N)
   integer IND(N)
@@ -222,11 +223,11 @@
    IND(j)=j
   enddo
 
-  if (n == 1) return
+  if(n == 1) return
   L=n/2+1
   ir=n
   100 continue
-   IF (l > 1) THEN
+   IF(l > 1) THEN
      l=l-1
      indx=ind(l)
      q=a(indx)
@@ -235,19 +236,19 @@
      q=a(indx)
      ind(ir)=ind(1)
      ir=ir-1
-     if (ir == 1) then
+     if(ir == 1) then
        ind(1)=indx
        return
      endif
    ENDIF
    i=l
    j=l+l
-  200    continue
-   IF (J <= IR) THEN
-      IF (J < IR) THEN
-         IF ( A(IND(j)) < A(IND(j+1)) ) j=j+1
+  200 continue
+   IF(J <= IR) THEN
+      IF(J < IR) THEN
+         IF(A(IND(j)) < A(IND(j+1))) j=j+1
       ENDIF
-      IF (q < A(IND(j))) THEN
+      IF(q < A(IND(j))) THEN
          IND(I)=IND(J)
          I=J
          J=J+J
@@ -275,12 +276,10 @@
 
   integer j
 
-  do J=1,N
-    W(j)=A(j)
-  enddo
+  W(:) = A(:)
 
   do J=1,N
-    A(j)=W(ind(j))
+    A(j) = W(ind(j))
   enddo
 
   end subroutine swap
@@ -298,12 +297,10 @@
 
   integer j
 
-  do J=1,N
-    W(j)=A(j)
-  enddo
+  W(:) = A(:)
 
   do J=1,N
-    A(j)=W(ind(j))
+    A(j) = W(ind(j))
   enddo
 
   end subroutine iswap

Modified: seismo/2D/SPECFEM2D/trunk/SPECFEM90/createnum_slow.f90
===================================================================
--- seismo/2D/SPECFEM2D/trunk/SPECFEM90/createnum_slow.f90	2004-12-14 15:22:01 UTC (rev 8432)
+++ seismo/2D/SPECFEM2D/trunk/SPECFEM90/createnum_slow.f90	2007-12-07 23:45:56 UTC (rev 8433)
@@ -1,13 +1,13 @@
 
 !========================================================================
 !
-!                   S P E C F E M 2 D  Version 5.0
+!                   S P E C F E M 2 D  Version 5.1
 !                   ------------------------------
 !
 !                         Dimitri Komatitsch
 !          Universite de Pau et des Pays de l'Adour, France
 !
-!                          (c) May 2004
+!                          (c) December 2004
 !
 !========================================================================
 
@@ -121,9 +121,9 @@
           else if(ngnodother == 4) then
             i2 = 1
             j2 = NGLLZ
-            else
-                  stop 'bad corner'
-            endif
+          else
+            stop 'bad corner'
+          endif
 
 ! affecter le meme numero
           ibool(i,j,numelem) = ibool(i2,j2,num2)
@@ -249,7 +249,7 @@
   enddo
 
 ! sortir de la recherche
-          goto 135
+        goto 135
 
       endif
     enddo

Modified: seismo/2D/SPECFEM2D/trunk/SPECFEM90/datim.f90
===================================================================
--- seismo/2D/SPECFEM2D/trunk/SPECFEM90/datim.f90	2004-12-14 15:22:01 UTC (rev 8432)
+++ seismo/2D/SPECFEM2D/trunk/SPECFEM90/datim.f90	2007-12-07 23:45:56 UTC (rev 8433)
@@ -1,13 +1,13 @@
 
 !========================================================================
 !
-!                   S P E C F E M 2 D  Version 5.0
+!                   S P E C F E M 2 D  Version 5.1
 !                   ------------------------------
 !
 !                         Dimitri Komatitsch
 !          Universite de Pau et des Pays de l'Adour, France
 !
-!                          (c) May 2004
+!                          (c) December 2004
 !
 !========================================================================
 

Modified: seismo/2D/SPECFEM2D/trunk/SPECFEM90/defarrays.f90
===================================================================
--- seismo/2D/SPECFEM2D/trunk/SPECFEM90/defarrays.f90	2004-12-14 15:22:01 UTC (rev 8432)
+++ seismo/2D/SPECFEM2D/trunk/SPECFEM90/defarrays.f90	2007-12-07 23:45:56 UTC (rev 8433)
@@ -1,13 +1,13 @@
 
 !========================================================================
 !
-!                   S P E C F E M 2 D  Version 5.0
+!                   S P E C F E M 2 D  Version 5.1
 !                   ------------------------------
 !
 !                         Dimitri Komatitsch
 !          Universite de Pau et des Pays de l'Adour, France
 !
-!                          (c) May 2004
+!                          (c) December 2004
 !
 !========================================================================
 

Modified: seismo/2D/SPECFEM2D/trunk/SPECFEM90/define_derivative_matrices.f90
===================================================================
--- seismo/2D/SPECFEM2D/trunk/SPECFEM90/define_derivative_matrices.f90	2004-12-14 15:22:01 UTC (rev 8432)
+++ seismo/2D/SPECFEM2D/trunk/SPECFEM90/define_derivative_matrices.f90	2007-12-07 23:45:56 UTC (rev 8433)
@@ -1,13 +1,13 @@
 
 !========================================================================
 !
-!                   S P E C F E M 2 D  Version 5.0
+!                   S P E C F E M 2 D  Version 5.1
 !                   ------------------------------
 !
 !                         Dimitri Komatitsch
 !          Universite de Pau et des Pays de l'Adour, France
 !
-!                          (c) May 2004
+!                          (c) December 2004
 !
 !========================================================================
 
@@ -34,7 +34,6 @@
 
   integer i1,i2,k1,k2
 
-! $$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$
 ! set up coordinates of the Gauss-Lobatto-Legendre points
   call zwgljd(xigll,wxgll,NGLLX,GAUSSALPHA,GAUSSBETA)
   call zwgljd(zigll,wzgll,NGLLZ,GAUSSALPHA,GAUSSBETA)

Modified: seismo/2D/SPECFEM2D/trunk/SPECFEM90/gll_library.f90
===================================================================
--- seismo/2D/SPECFEM2D/trunk/SPECFEM90/gll_library.f90	2004-12-14 15:22:01 UTC (rev 8432)
+++ seismo/2D/SPECFEM2D/trunk/SPECFEM90/gll_library.f90	2007-12-07 23:45:56 UTC (rev 8433)
@@ -20,14 +20,14 @@
   integer i
 
   f3 = zero
-  apb   = alpha+beta
-  if (n == 0) then
-   endw1 = zero
-   return
+  apb = alpha+beta
+  if(n == 0) then
+    endw1 = zero
+    return
   endif
-  f1   = gammaf(alpha+two)*gammaf(beta+one)/gammaf(apb+three)
-  f1   = f1*(apb+two)*two**(apb+two)/two
-  if (n == 1) then
+  f1 = gammaf(alpha+two)*gammaf(beta+one)/gammaf(apb+three)
+  f1 = f1*(apb+two)*two**(apb+two)/two
+  if(n == 1) then
    endw1 = f1
    return
   endif
@@ -36,7 +36,7 @@
   fint2 = gammaf(alpha+two)*gammaf(beta+two)/gammaf(apb+four)
   fint2 = fint2*two**(apb+three)
   f2    = (-two*(beta+two)*fint1 + (apb+four)*fint2) * (apb+three)/four
-  if (n == 2) then
+  if(n == 2) then
    endw1 = f2
    return
   endif

Modified: seismo/2D/SPECFEM2D/trunk/SPECFEM90/gmat01.f90
===================================================================
--- seismo/2D/SPECFEM2D/trunk/SPECFEM90/gmat01.f90	2004-12-14 15:22:01 UTC (rev 8432)
+++ seismo/2D/SPECFEM2D/trunk/SPECFEM90/gmat01.f90	2007-12-07 23:45:56 UTC (rev 8433)
@@ -1,13 +1,13 @@
 
 !========================================================================
 !
-!                   S P E C F E M 2 D  Version 5.0
+!                   S P E C F E M 2 D  Version 5.1
 !                   ------------------------------
 !
 !                         Dimitri Komatitsch
 !          Universite de Pau et des Pays de l'Adour, France
 !
-!                          (c) May 2004
+!                          (c) December 2004
 !
 !========================================================================
 
@@ -59,29 +59,13 @@
 ! Poisson's ratio must be between -1 and +1/2
       if (poisson < -1.d0 .or. poisson > 0.5d0) stop 'Poisson''s ratio out of range'
 
-!---- materiau isotrope, module de Young et coefficient de Poisson donnes
-   else if(indic == 1) then
-      young = val1
-      poisson = val2
-! Poisson's ratio must be between -1 and +1/2
-      if (poisson < -1.d0 .or. poisson > 0.5d0) stop 'Poisson''s ratio out of range'
-      a2mu  = young/(one+poisson)
-      amu   = half*a2mu
-      alam  = a2mu*poisson/(one-two*poisson)
-      Kmod  = alam + a2mu
-      Kvol  = alam + a2mu/3.d0
-      cp    = dsqrt((Kvol + 4.d0*amu/3.d0)/denst)
-      cs    = dsqrt(amu/denst)
-
 !---- materiau anisotrope, c11, c13, c33 et c44 donnes en Pascal
-   else if(indic == 2) then
+   else
       c11 = val1
       c13 = val2
       c33 = val3
       c44 = val4
 
-   else
-      stop 'Improper value while reading material sets'
    endif
 
 !

Modified: seismo/2D/SPECFEM2D/trunk/SPECFEM90/lagrange_poly.f90
===================================================================
--- seismo/2D/SPECFEM2D/trunk/SPECFEM90/lagrange_poly.f90	2004-12-14 15:22:01 UTC (rev 8432)
+++ seismo/2D/SPECFEM2D/trunk/SPECFEM90/lagrange_poly.f90	2007-12-07 23:45:56 UTC (rev 8433)
@@ -1,21 +1,22 @@
 
 !========================================================================
 !
-!                   S P E C F E M 2 D  Version 5.0
+!                   S P E C F E M 2 D  Version 5.1
 !                   ------------------------------
 !
 !                         Dimitri Komatitsch
 !          Universite de Pau et des Pays de l'Adour, France
 !
-!                          (c) May 2004
+!                          (c) December 2004
 !
 !========================================================================
 
   double precision function hgll(I,Z,ZGLL,NZ)
+
 !-------------------------------------------------------------
 !
-!     Compute the value of the Lagrangian interpolant L through
-!     the NZ Gauss-Lobatto Legendre points ZGLL at the point Z
+!  Compute the value of the Lagrangian interpolant L through
+!  the NZ Gauss-Lobatto Legendre points ZGLL at point Z
 !
 !-------------------------------------------------------------
 

Modified: seismo/2D/SPECFEM2D/trunk/SPECFEM90/plotgll.f90
===================================================================
--- seismo/2D/SPECFEM2D/trunk/SPECFEM90/plotgll.f90	2004-12-14 15:22:01 UTC (rev 8432)
+++ seismo/2D/SPECFEM2D/trunk/SPECFEM90/plotgll.f90	2007-12-07 23:45:56 UTC (rev 8433)
@@ -1,13 +1,13 @@
 
 !========================================================================
 !
-!                   S P E C F E M 2 D  Version 5.0
+!                   S P E C F E M 2 D  Version 5.1
 !                   ------------------------------
 !
 !                         Dimitri Komatitsch
 !          Universite de Pau et des Pays de l'Adour, France
 !
-!                          (c) May 2004
+!                          (c) December 2004
 !
 !========================================================================
 
@@ -39,7 +39,7 @@
   write(iout,*) 'Generating gnuplot meshes...'
   write(iout,*)
 
-! create non empty files for the case of 4-nodes elements
+! create non empty files for the case of 4-node elements
 
   name='macros1.gnu'
   open(unit=30,file=name,status='unknown')
@@ -125,7 +125,7 @@
 
   if(ngnod == 4) then
 !
-!----  4-noded rectangular element
+!----  4-node rectangular element
 !
 
 ! draw the edges of the element using one color
@@ -145,7 +145,7 @@
   else
 
 !
-!----  9-noded rectangular element
+!----  9-node rectangular element
 !
 
 ! draw the edges of the element using one color
@@ -223,10 +223,6 @@
   write(20,*) 'pause -1 "Hit any key to exit..."'
   close(20)
 
-!
-!----
-!
-
  10 format('')
  15 format(e10.5,1x,e10.5)
 

Modified: seismo/2D/SPECFEM2D/trunk/SPECFEM90/plotpost.f90
===================================================================
--- seismo/2D/SPECFEM2D/trunk/SPECFEM90/plotpost.f90	2004-12-14 15:22:01 UTC (rev 8432)
+++ seismo/2D/SPECFEM2D/trunk/SPECFEM90/plotpost.f90	2007-12-07 23:45:56 UTC (rev 8433)
@@ -1,13 +1,13 @@
 
 !========================================================================
 !
-!                   S P E C F E M 2 D  Version 5.0
+!                   S P E C F E M 2 D  Version 5.1
 !                   ------------------------------
 !
 !                         Dimitri Komatitsch
 !          Universite de Pau et des Pays de l'Adour, France
 !
-!                          (c) May 2004
+!                          (c) December 2004
 !
 !========================================================================
 

Modified: seismo/2D/SPECFEM2D/trunk/SPECFEM90/positrec.f90
===================================================================
--- seismo/2D/SPECFEM2D/trunk/SPECFEM90/positrec.f90	2004-12-14 15:22:01 UTC (rev 8432)
+++ seismo/2D/SPECFEM2D/trunk/SPECFEM90/positrec.f90	2007-12-07 23:45:56 UTC (rev 8433)
@@ -1,13 +1,13 @@
 
 !========================================================================
 !
-!                   S P E C F E M 2 D  Version 5.0
+!                   S P E C F E M 2 D  Version 5.1
 !                   ------------------------------
 !
 !                         Dimitri Komatitsch
 !          Universite de Pau et des Pays de l'Adour, France
 !
-!                          (c) May 2004
+!                          (c) December 2004
 !
 !========================================================================
 

Modified: seismo/2D/SPECFEM2D/trunk/SPECFEM90/positsource.f90
===================================================================
--- seismo/2D/SPECFEM2D/trunk/SPECFEM90/positsource.f90	2004-12-14 15:22:01 UTC (rev 8432)
+++ seismo/2D/SPECFEM2D/trunk/SPECFEM90/positsource.f90	2007-12-07 23:45:56 UTC (rev 8433)
@@ -1,13 +1,13 @@
 
 !========================================================================
 !
-!                   S P E C F E M 2 D  Version 5.0
+!                   S P E C F E M 2 D  Version 5.1
 !                   ------------------------------
 !
 !                         Dimitri Komatitsch
 !          Universite de Pau et des Pays de l'Adour, France
 !
-!                          (c) May 2004
+!                          (c) December 2004
 !
 !========================================================================
 

Modified: seismo/2D/SPECFEM2D/trunk/SPECFEM90/q49shape.f90
===================================================================
--- seismo/2D/SPECFEM2D/trunk/SPECFEM90/q49shape.f90	2004-12-14 15:22:01 UTC (rev 8432)
+++ seismo/2D/SPECFEM2D/trunk/SPECFEM90/q49shape.f90	2007-12-07 23:45:56 UTC (rev 8433)
@@ -1,13 +1,13 @@
 
 !========================================================================
 !
-!                   S P E C F E M 2 D  Version 5.0
+!                   S P E C F E M 2 D  Version 5.1
 !                   ------------------------------
 !
 !                         Dimitri Komatitsch
 !          Universite de Pau et des Pays de l'Adour, France
 !
-!                          (c) May 2004
+!                          (c) December 2004
 !
 !========================================================================
 
@@ -54,7 +54,7 @@
 !
   if(ngnod == 4) then
 !
-!----    4-noded rectangular element
+!----    4-node rectangular element
 !
  do l2 = 1,NGLLZ
 
@@ -92,7 +92,7 @@
 
   else if(ngnod == 9) then
 !
-!----    9-noded rectangular element
+!----    9-node rectangular element
 !
  do l2 = 1,NGLLZ
 

Modified: seismo/2D/SPECFEM2D/trunk/SPECFEM90/q49spec.f90
===================================================================
--- seismo/2D/SPECFEM2D/trunk/SPECFEM90/q49spec.f90	2004-12-14 15:22:01 UTC (rev 8432)
+++ seismo/2D/SPECFEM2D/trunk/SPECFEM90/q49spec.f90	2007-12-07 23:45:56 UTC (rev 8433)
@@ -1,13 +1,13 @@
 
 !========================================================================
 !
-!                   S P E C F E M 2 D  Version 5.0
+!                   S P E C F E M 2 D  Version 5.1
 !                   ------------------------------
 !
 !                         Dimitri Komatitsch
 !          Universite de Pau et des Pays de l'Adour, France
 !
-!                          (c) May 2004
+!                          (c) December 2004
 !
 !========================================================================
 
@@ -126,7 +126,7 @@
 !
   if(ngnod == 4) then
 !
-!----    4-noded rectangular element
+!----    4-node rectangular element
 !
  do l2 = 1,iptsdisp
 
@@ -154,7 +154,7 @@
 
   else if(ngnod == 9) then
 !
-!----    9-noded rectangular element
+!----    9-node rectangular element
 !
  do l2 = 1,iptsdisp
 

Modified: seismo/2D/SPECFEM2D/trunk/SPECFEM90/specfem2D.f90
===================================================================
--- seismo/2D/SPECFEM2D/trunk/SPECFEM90/specfem2D.f90	2004-12-14 15:22:01 UTC (rev 8432)
+++ seismo/2D/SPECFEM2D/trunk/SPECFEM90/specfem2D.f90	2007-12-07 23:45:56 UTC (rev 8433)
@@ -1,30 +1,39 @@
 
 !========================================================================
 !
-!                   S P E C F E M 2 D  Version 5.0
+!                   S P E C F E M 2 D  Version 5.1
 !                   ------------------------------
 !
 !                         Dimitri Komatitsch
 !          Universite de Pau et des Pays de l'Adour, France
 !
-!                          (c) May 2004
+!                          (c) December 2004
 !
 !========================================================================
 
-!========================================================================
+!====================================================================================
 !
-!  An explicit 2D spectral element solver for the elastic wave equation
+! An explicit 2D spectral element solver for the anelastic anisotropic wave equation
 !
-!========================================================================
+!====================================================================================
 
-! version 5.0 : - got rid of useless routines, suppressed commons etc.
+!
+! version 5.1, December 2004 : more general mesher with any number of curved layers
+!
+! version 5.0, May 2004 :
+!               - got rid of useless routines, suppressed commons etc.
 !               - weak formulation based explicitly on stress tensor
 !               - implementation of full anisotropy
 !               - implementation of attenuation based on memory variables
-
-! version 5.0 is based on SPECFEM2D version 4.2
-! (c) June 1998 by Dimitri Komatitsch, Harvard University, USA
+!
+! based on SPECFEM2D version 4.2, June 1998
+! (c) by Dimitri Komatitsch, Harvard University, USA
 ! and Jean-Pierre Vilotte, Institut de Physique du Globe de Paris, France
+!
+! itself based on SPECFEM2D version 1.0, 1995
+! (c) by Dimitri Komatitsch and Jean-Pierre Vilotte,
+! Institut de Physique du Globe de Paris, France
+!
 
   program specfem2D
 
@@ -45,7 +54,7 @@
   logical anyabs
 
   integer i,j,it,irec,iglobrec,ipoin,ip,id,ip1,ip2,in,nnum
-  integer nbpoin,inump,n,npoinext,netyp,ispec,npoin,npgeo,iglob
+  integer nbpoin,inump,n,npoinext,ispec,npoin,npgeo,iglob
 
   double precision valux,valuz,rhoextread,vpextread,vsextread
   double precision cpl,csl,rhol
@@ -102,7 +111,7 @@
   integer, dimension(:,:), allocatable  :: knods
   integer, dimension(:), allocatable :: kmato,numabs
 
-  integer ie,k,isource,iexplo
+  integer ie,k,iexplo
 
   integer ielems,iglobsource
   double precision f0,t0,factor,a,angleforce,ricker,displnorm_all
@@ -237,14 +246,12 @@
 !----  read source information
 !
   read(IIN,40) datlin
-  read(IIN,*) (gltfu(k), k=1,9)
+  read(IIN,*) (gltfu(k), k=2,9)
 
 !
 !-----  check the input
 !
  if(.not. initialfield) then
-   if(nint(gltfu(1)) /= 6) stop 'Wrong function number in getltf !'
-   isource = nint(gltfu(1))
    iexplo = nint(gltfu(2))
    if (iexplo == 1) then
      write(IOUT,212) (gltfu(k), k=3,8)
@@ -255,12 +262,8 @@
    endif
  endif
 
-!
 !-----  convert angle from degrees to radians
-!
-   isource = nint(gltfu(1))
-   iexplo = nint(gltfu(2))
-   if(isource >= 4 .and. isource <= 6 .and. iexplo == 1) gltfu(8) = gltfu(8) * pi / 180.d0
+  gltfu(8) = gltfu(8) * pi / 180.d0
 
 !
 !---- read receiver locations
@@ -292,7 +295,7 @@
 !---- read the basic properties of the spectral elements
 !
   read(IIN ,40) datlin
-  read(IIN ,*) netyp,numat,ngnod,nspec,iptsdisp,nelemabs
+  read(IIN ,*) numat,ngnod,nspec,iptsdisp,nelemabs
 
 !
 !---- allocate arrays

Modified: seismo/2D/SPECFEM2D/trunk/SPECFEM90/write_seismograms.f90
===================================================================
--- seismo/2D/SPECFEM2D/trunk/SPECFEM90/write_seismograms.f90	2004-12-14 15:22:01 UTC (rev 8432)
+++ seismo/2D/SPECFEM2D/trunk/SPECFEM90/write_seismograms.f90	2007-12-07 23:45:56 UTC (rev 8433)
@@ -1,13 +1,13 @@
 
 !========================================================================
 !
-!                   S P E C F E M 2 D  Version 5.0
+!                   S P E C F E M 2 D  Version 5.1
 !                   ------------------------------
 !
 !                         Dimitri Komatitsch
 !          Universite de Pau et des Pays de l'Adour, France
 !
-!                          (c) May 2004
+!                          (c) December 2004
 !
 !========================================================================
 



More information about the cig-commits mailing list