[cig-commits] r22662 - in seismo/3D/SPECFEM3D_GLOBE: branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D trunk/src/meshfem3D

dkomati1 at geodynamics.org dkomati1 at geodynamics.org
Wed Jul 24 08:33:12 PDT 2013


Author: dkomati1
Date: 2013-07-24 08:33:12 -0700 (Wed, 24 Jul 2013)
New Revision: 22662

Modified:
   seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/model_1dref.f90
   seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/model_aniso_mantle.f90
   seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/model_attenuation.f90
   seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/model_crust.f90
   seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/model_crustmaps.f90
   seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/model_epcrust.f90
   seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/model_eucrust.f90
   seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/model_gapp2.f90
   seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/model_gll.f90
   seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/model_heterogen_mantle.f90
   seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/model_jp3d.f90
   seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/model_ppm.f90
   seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/model_s20rts.f90
   seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/model_s362ani.f90
   seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/model_s40rts.f90
   seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/rules.mk
   seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/save_arrays_solver.f90
   seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/model_1066a.f90
   seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/model_1dref.f90
   seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/model_ak135.f90
   seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/model_aniso_mantle.f90
   seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/model_atten3D_QRFSI12.f90
   seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/model_attenuation.f90
   seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/model_crust.f90
   seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/model_crustmaps.f90
   seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/model_epcrust.f90
   seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/model_heterogen_mantle.f90
   seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/model_jp3d.f90
   seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/model_s20rts.f90
Log:
done with all the easy merges in src/meshfem3D


Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/model_1dref.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/model_1dref.f90	2013-07-24 14:30:04 UTC (rev 22661)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/model_1dref.f90	2013-07-24 15:33:12 UTC (rev 22662)
@@ -125,7 +125,6 @@
 ! shear quality factor Qmu
 ! bulk quality factor Qkappa
 
-
   double precision :: x,rho,vpv,vph,vsv,vsh,eta,Qmu,Qkappa
   integer :: iregion_code
   logical :: CRUSTAL

Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/model_aniso_mantle.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/model_aniso_mantle.f90	2013-07-24 14:30:04 UTC (rev 22661)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/model_aniso_mantle.f90	2013-07-24 15:33:12 UTC (rev 22662)
@@ -58,12 +58,11 @@
 ! standard routine to setup model
 
   use model_aniso_mantle_par
+  use mpi
 
   implicit none
 
   include "constants.h"
-  ! standard include of the MPI library
-  include 'mpif.h'
 
   integer :: myrank
 

Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/model_attenuation.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/model_attenuation.f90	2013-07-24 14:30:04 UTC (rev 22661)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/model_attenuation.f90	2013-07-24 15:33:12 UTC (rev 22662)
@@ -65,11 +65,11 @@
 
 ! standard routine to setup model
 
+  use mpi
+
   implicit none
 
   include "constants.h"
-  ! standard include of the MPI library
-  include 'mpif.h'
 
 ! model_attenuation_variables
   type model_attenuation_variables
@@ -537,6 +537,8 @@
 
   implicit none
 
+  include "constants.h"
+
   integer n
   double precision tau_s(n)
   double precision min_period, max_period
@@ -544,7 +546,6 @@
   double precision exp1, exp2
   double precision dexp
   integer i
-  double precision, parameter :: PI = 3.14159265358979d0
 
   f1 = 1.0d0 / max_period
   f2 = 1.0d0 / min_period
@@ -567,6 +568,8 @@
 
   implicit none
 
+  include "constants.h"
+
 ! attenuation_simplex_variables
   type attenuation_simplex_variables
     sequence
@@ -596,7 +599,6 @@
   integer i, iterations, err,prnt
   double precision f1, f2, exp1,exp2,dexp, min_value
   double precision, allocatable, dimension(:) :: f
-  double precision, parameter :: PI = 3.14159265358979d0
   integer, parameter :: nf = 100
   double precision, external :: attenuation_eval
 
@@ -768,6 +770,8 @@
 
   implicit none
 
+  include "constants.h"
+
   ! Input
   integer nf, nsls
   double precision, dimension(nf)   :: f
@@ -776,21 +780,17 @@
   double precision, dimension(nf)   :: A,B
 
   integer i,j
-  double precision w, pi, demon
+  double precision w, demon
 
-  PI = 3.14159265358979d0
-
   A(:) = 1.0d0 -  nsls*1.0d0
   B(:) = 0.0d0
   do i = 1,nf
      w = 2.0d0 * PI * 10**f(i)
      do j = 1,nsls
-!        write(*,*)j,tau_s(j),tau_e(j)
         demon = 1.0d0 + w**2 * tau_s(j)**2
         A(i) = A(i) + ((1.0d0 + (w**2 * tau_e(j) * tau_s(j)))/ demon)
         B(i) = B(i) + ((w * (tau_e(j) - tau_s(j))) / demon)
      enddo
-!     write(*,*)A(i),B(i),10**f(i)
   enddo
 
   end subroutine attenuation_maxwell
@@ -1260,209 +1260,3 @@
 
   end subroutine qsort_local
 
-
-!
-!-------------------------------------------------------------------------------------------------
-!
-! unused routines...
-!
-!
-!  subroutine model_attenuation_1D_PREM(x, Qmu)
-!
-!! x is the radius from 0 to 1 where 0 is the center and 1 is the surface
-!! This version is for 1D PREM.
-!
-!  implicit none
-!
-!  include 'constants.h'
-!!  integer iflag
-!  double precision r, x, Qmu,RICB,RCMB, &
-!      RTOPDDOUBLEPRIME,R600,R670,R220,R771,R400,R80, ROCEAN, RMOHO, RMIDDLE_CRUST
-!  double precision Qkappa
-!
-!  r = x * R_EARTH
-!
-!  ROCEAN = 6368000.d0
-!  RMIDDLE_CRUST = 6356000.d0
-!  RMOHO = 6346600.d0
-!  R80  = 6291000.d0
-!  R220 = 6151000.d0
-!  R400 = 5971000.d0
-!  R600 = 5771000.d0
-!  R670 = 5701000.d0
-!  R771 = 5600000.d0
-!  RTOPDDOUBLEPRIME = 3630000.d0
-!  RCMB = 3480000.d0
-!  RICB = 1221000.d0
-!
-!! PREM
-!!
-!!--- inner core
-!!
-!  if(r >= 0.d0 .and. r <= RICB) then
-!     Qmu=84.6d0
-!     Qkappa=1327.7d0
-!!
-!!--- outer core
-!!
-!  else if(r > RICB .and. r <= RCMB) then
-!     Qmu=0.0d0
-!     Qkappa=57827.0d0
-!     if(RCMB - r < r - RICB) then
-!        Qmu = 312.0d0  ! CMB
-!     else
-!        Qmu = 84.6d0   ! ICB
-!     endif
-!!
-!!--- D" at the base of the mantle
-!!
-!  else if(r > RCMB .and. r <= RTOPDDOUBLEPRIME) then
-!     Qmu=312.0d0
-!     Qkappa=57827.0d0
-!!
-!!--- mantle: from top of D" to d670
-!!
-!  else if(r > RTOPDDOUBLEPRIME .and. r <= R771) then
-!     Qmu=312.0d0
-!     Qkappa=57827.0d0
-!  else if(r > R771 .and. r <= R670) then
-!     Qmu=312.0d0
-!     Qkappa=57827.0d0
-!!
-!!--- mantle: above d670
-!!
-!  else if(r > R670 .and. r <= R600) then
-!     Qmu=143.0d0
-!     Qkappa=57827.0d0
-!  else if(r > R600 .and. r <= R400) then
-!     Qmu=143.0d0
-!     Qkappa=57827.0d0
-!  else if(r > R400 .and. r <= R220) then
-!     Qmu=143.0d0
-!     Qkappa=57827.0d0
-!  else if(r > R220 .and. r <= R80) then
-!     Qmu=80.0d0
-!     Qkappa=57827.0d0
-!  else if(r > R80) then
-!     Qmu=600.0d0
-!     Qkappa=57827.0d0
-!  endif
-!
-!! Since R80 may be changed, we use radius to decide the attenuation region
-!! rather than doubling flag
-!
-!  ! We determine the attenuation value here dependent on the doubling flag and
-!  ! which region we are sitting in. The radius reported is not accurate for
-!  ! determination of which region we are actually in, whereas the idoubling flag is
-!!  if(iflag == IFLAG_INNER_CORE_NORMAL .or. iflag == IFLAG_MIDDLE_CENTRAL_CUBE .or. &
-!!       iflag == IFLAG_BOTTOM_CENTRAL_CUBE .or. iflag == IFLAG_TOP_CENTRAL_CUBE .or. &
-!!       iflag == IFLAG_IN_FICTITIOUS_CUBE) then
-!!     Qmu =  84.6d0
-!!     Qkappa = 1327.7d0
-!!  else if(iflag == IFLAG_OUTER_CORE_NORMAL) then
-!!     Qmu = 0.0d0
-!!     Qkappa = 57827.0d0
-!!  else if(iflag == IFLAG_MANTLE_NORMAL) then ! D'' to 670 km
-!!     Qmu = 312.0d0
-!!     Qkappa = 57827.0d0
-!!  else if(iflag == IFLAG_670_220) then
-!!     Qmu=143.0d0
-!!     Qkappa = 57827.0d0
-!!  else if(iflag == IFLAG_220_80) then
-!!     Qmu=80.0d0
-!!     Qkappa = 57827.0d0
-!!  else if(iflag == IFLAG_80_MOHO) then
-!!     Qmu=600.0d0
-!!     Qkappa = 57827.0d0
-!!  else if(iflag == IFLAG_CRUST) then
-!!     Qmu=600.0d0
-!!     Qkappa = 57827.0d0
-!!  else
-!!     write(*,*)'iflag:',iflag
-!!     call exit_MPI_without_rank('Invalid idoubling flag in attenuation_model_1D_prem from get_model()')
-!!  endif
-!
-!  end subroutine model_attenuation_1D_PREM
-!
-!!
-!!-------------------------------------------------------------------------------------------------
-!!
-!
-!! get 1D REF attenuation model according to radius
-!  subroutine model_attenuation_1D_REF(x, Qmu)
-!
-!! x in the radius from 0 to 1 where 0 is the center and 1 is the surface
-!! This version is for 1D REF.
-!
-!  implicit none
-!
-!  include 'constants.h'
-!
-!  double precision r, x, Qmu,RICB,RCMB, &
-!      RTOPDDOUBLEPRIME,R600,R670,R220,R771,R400,R80, ROCEAN, RMOHO, RMIDDLE_CRUST
-!  double precision Qkappa
-!
-!  r = x * R_EARTH
-!
-!  ROCEAN = 6368000.d0
-!  RMIDDLE_CRUST = 6356000.d0
-!  RMOHO = 6346600.d0
-!  R80  = 6291000.d0
-!  R220 = 6151000.d0
-!  R400 = 5961000.d0
-!  R600 = 5771000.d0
-!  R670 = 5721000.d0
-!  R771 = 5600000.d0
-!  RTOPDDOUBLEPRIME = 3630000.d0
-!  RCMB = 3479958.d0
-!  RICB = 1221491.d0
-!
-!! REF model
-!!
-!!--- inner core
-!!
-!  if(r >= 0.d0 .and. r <= RICB) then
-!     Qmu=104.0d0
-!     Qkappa=1327.6d0
-!
-!!--- outer core
-!!
-!  else if(r > RICB .and. r <= RCMB) then
-!     Qmu=0.0d0
-!     Qkappa=57822.5d0
-!     if(RCMB - r < r - RICB) then
-!        Qmu = 355.0d0  ! CMB
-!     else
-!        Qmu = 104.0d0   ! ICB
-!     endif
-!
-!!--- D" at the base of the mantle
-!!
-!  else if(r > RCMB .and. r <= RTOPDDOUBLEPRIME) then
-!     Qmu=355.0d0
-!     Qkappa=57822.5d0
-!
-!!--- mantle: from top of D" to d670
-!!
-!  else if(r > RTOPDDOUBLEPRIME .and. r <= R670) then
-!     Qmu=355.0d0
-!     Qkappa=57822.5d0
-!
-!!--- mantle: above d670
-!!
-!  else if(r > R670 .and. r <= R220) then
-!     Qmu=165.0d0
-!     Qkappa=943.0d0
-!  else if(r > R220 .and. r <= R80) then
-!     Qmu=70.0d0
-!     Qkappa=943.0d0
-!  else if(r > R80.and. r<=RMOHO) then
-!     Qmu=191.0d0
-!     Qkappa=943.0d0
-!  else if (r > RMOHO) then
-!     Qmu=300.0d0
-!     Qkappa=57822.5d0
-!  endif
-!
-!  end subroutine model_attenuation_1D_REF
-!

Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/model_crust.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/model_crust.f90	2013-07-24 14:30:04 UTC (rev 22661)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/model_crust.f90	2013-07-24 15:33:12 UTC (rev 22662)
@@ -58,12 +58,12 @@
 
 ! standard routine to setup model
 
+  use mpi
   use model_crust_par
 
   implicit none
 
   include "constants.h"
-  include 'mpif.h'
 
   integer :: myrank
   integer :: ier

Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/model_crustmaps.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/model_crustmaps.f90	2013-07-24 14:30:04 UTC (rev 22661)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/model_crustmaps.f90	2013-07-24 15:33:12 UTC (rev 22662)
@@ -62,12 +62,12 @@
 
 ! standard routine to setup model
 
+  use mpi
   use model_crustmaps_par
 
   implicit none
 
   include "constants.h"
-  include 'mpif.h'
 
   integer :: myrank
 

Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/model_epcrust.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/model_epcrust.f90	2013-07-24 14:30:04 UTC (rev 22661)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/model_epcrust.f90	2013-07-24 15:33:12 UTC (rev 22662)
@@ -61,12 +61,12 @@
 
   subroutine model_epcrust_broadcast(myrank)
 
+  use mpi
   use model_epcrust_par
 
   implicit none
 
   include "constants.h"
-  include 'mpif.h'
 
   integer:: myrank,ier
 

Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/model_eucrust.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/model_eucrust.f90	2013-07-24 14:30:04 UTC (rev 22661)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/model_eucrust.f90	2013-07-24 15:33:12 UTC (rev 22662)
@@ -51,12 +51,12 @@
 
 ! standard routine to setup model
 
+  use mpi
   use model_eucrust_par
 
   implicit none
 
   include "constants.h"
-  include 'mpif.h'
 
   integer :: myrank
   integer :: ier

Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/model_gapp2.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/model_gapp2.f90	2013-07-24 14:30:04 UTC (rev 22661)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/model_gapp2.f90	2013-07-24 15:33:12 UTC (rev 22662)
@@ -53,11 +53,12 @@
 ! standard routine to setup model
 
   use gapp2_mantle_model_constants
+  use mpi
 
   implicit none
+
   include "constants.h"
-  ! standard include of the MPI library
-  include 'mpif.h'
+
   integer :: myrank
   integer :: ier
 

Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/model_gll.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/model_gll.f90	2013-07-24 14:30:04 UTC (rev 22661)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/model_gll.f90	2013-07-24 15:33:12 UTC (rev 22662)
@@ -38,12 +38,11 @@
 ! standard routine to setup model
 
   use meshfem3D_models_par,only: TRANSVERSE_ISOTROPY
+  use mpi
 
   implicit none
 
   include "constants.h"
-  ! standard include of the MPI library
-  include 'mpif.h'
   include "precision.h"
 
   ! GLL model_variables

Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/model_heterogen_mantle.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/model_heterogen_mantle.f90	2013-07-24 14:30:04 UTC (rev 22661)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/model_heterogen_mantle.f90	2013-07-24 15:33:12 UTC (rev 22662)
@@ -49,13 +49,12 @@
 
 ! standard routine to setup model
 
+  use mpi
   use model_heterogen_mantle_par
 
   implicit none
 
   include "constants.h"
-  ! standard include of the MPI library
-  include 'mpif.h'
 
   integer :: myrank
 

Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/model_jp3d.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/model_jp3d.f90	2013-07-24 14:30:04 UTC (rev 22661)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/model_jp3d.f90	2013-07-24 15:33:12 UTC (rev 22662)
@@ -125,13 +125,12 @@
 
 ! standard routine to setup model
 
+  use mpi
   use model_jp3d_par
 
   implicit none
 
   include "constants.h"
-  ! standard include of the MPI library
-  include 'mpif.h'
 
   integer :: myrank
   integer :: ier

Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/model_ppm.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/model_ppm.f90	2013-07-24 14:30:04 UTC (rev 22661)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/model_ppm.f90	2013-07-24 15:33:12 UTC (rev 22662)
@@ -96,13 +96,12 @@
 
 ! standard routine to setup model
 
+  use mpi
   use model_ppm_par
 
   implicit none
 
   include "constants.h"
-  ! standard include of the MPI library
-  include 'mpif.h'
 
   integer :: myrank
   integer :: ier

Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/model_s20rts.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/model_s20rts.f90	2013-07-24 14:30:04 UTC (rev 22661)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/model_s20rts.f90	2013-07-24 15:33:12 UTC (rev 22662)
@@ -57,18 +57,16 @@
 !--------------------------------------------------------------------------------------------------
 !
 
-
   subroutine model_s20rts_broadcast(myrank)
 
 ! standard routine to setup model
 
+  use mpi
   use model_s20rts_par
 
   implicit none
 
   include "constants.h"
-  ! standard include of the MPI library
-  include 'mpif.h'
 
   integer :: myrank
 

Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/model_s362ani.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/model_s362ani.f90	2013-07-24 14:30:04 UTC (rev 22661)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/model_s362ani.f90	2013-07-24 15:33:12 UTC (rev 22662)
@@ -83,13 +83,12 @@
 
 ! standard routine to setup model
 
+  use mpi
   use model_s362ani_par
 
   implicit none
 
   include "constants.h"
-  ! standard include of the MPI library
-  include 'mpif.h'
 
   integer :: myrank
   integer :: THREE_D_MODEL
@@ -1077,13 +1076,15 @@
 !-------------------------------------------------------------------------------------------------
 !
 
-
   subroutine splcon(xlat,xlon,nver,verlat,verlon,verrad,ncon,icon,con)
 
   implicit none
 
-  integer :: ncon,nver
+  include "constants.h"
 
+  integer, intent(in) :: nver
+  integer, intent(out) :: ncon
+
 ! Daniel Peter: original define
 !
 !  real(kind=4) verlat(1)
@@ -1094,39 +1095,39 @@
 !  real(kind=4) con(1)
 
 ! Daniel Peter: avoiding out-of-bounds errors
-  real(kind=4) verlat(nver)
-  real(kind=4) verlon(nver)
-  real(kind=4) verrad(nver)
+  real(kind=4), intent(in) :: verlat(nver)
+  real(kind=4), intent(in) :: verlon(nver)
+  real(kind=4), intent(in) :: verrad(nver)
 
   integer, parameter :: maxver=1000
-  integer icon(maxver)
-  real(kind=4) con(maxver)
+  integer, intent(out) :: icon(maxver)
+  real(kind=4), intent(out) :: con(maxver)
 
   double precision dd
   double precision rn
   double precision dr
-  double precision xrad
   double precision ver8
   double precision xla8
 
   integer :: iver
   real(kind=4) :: xlat,xlon
 
-  xrad=3.14159265358979/180.d0
-
   ncon=0
 
   do iver=1,nver
     if(xlat > verlat(iver)-2.*verrad(iver)) then
       if(xlat < verlat(iver)+2.*verrad(iver)) then
-        ver8=xrad*(verlat(iver))
-        xla8=xrad*(xlat)
-        dd=sin(ver8)*sin(xla8)
-        dd=dd+cos(ver8)*cos(xla8)* cos(xrad*(xlon-verlon(iver)))
-        dd=acos(dd)/xrad
+        ver8=DEGREES_TO_RADIANS*(verlat(iver))
+        xla8=DEGREES_TO_RADIANS*(xlat)
+        dd=sin(ver8)*sin(xla8) + cos(ver8)*cos(xla8)* cos(DEGREES_TO_RADIANS*(xlon-verlon(iver)))
+        dd=acos(dd) * RADIANS_TO_DEGREES
         if(dd > (verrad(iver))*2.d0) then
         else
           ncon=ncon+1
+
+!! DK DK added this safety test
+          if(ncon > maxver) stop 'error: ncon > maxver in splcon()'
+
           icon(ncon)=iver
           rn=dd/(verrad(iver))
           dr=rn-1.d0
@@ -1201,15 +1202,21 @@
         call ylm(y,x,lmax,ylmcof(1,ihpa),wk1,wk2,wk3)
       else if(itypehpa(ihpa) == 2) then
         numcof=numcoe(ihpa)
+
 ! originally called
-!        call splcon(y,x,numcof,xlaspl(1,ihpa), &
-!              xlospl(1,ihpa),radspl(1,ihpa), &
-!              nconpt(ihpa),iconpt(1,ihpa),conpt(1,ihpa))
+         call splcon(y,x,numcof,xlaspl(1,ihpa), &
+               xlospl(1,ihpa),radspl(1,ihpa), &
+               nconpt(ihpa),iconpt(1,ihpa),conpt(1,ihpa))
 
 ! making sure of array bounds
-        call splcon(y,x,numcof,xlaspl(1:numcof,ihpa), &
-              xlospl(1:numcof,ihpa),radspl(1:numcof,ihpa), &
-              nconpt(ihpa),iconpt(1,ihpa),conpt(1,ihpa))
+!! note from DK DK, July 2013: the code below is correct but it is exactly the same as the code
+!! note from DK DK, July 2013: above because xlaspl(1:numcof,ihpa) is contiguous in memory, and thus
+!! note from DK DK, July 2013: using a pointer to it such as xlaspl(1,ihpa) in the call seems fine as well;
+!! note from DK DK, July 2013: and some compilers could create expensive and useless memory copies for each call
+!! note from DK DK, July 2013: with the new syntax below
+!       call splcon(y,x,numcof,xlaspl(1:numcof,ihpa), &
+!             xlospl(1:numcof,ihpa),radspl(1:numcof,ihpa), &
+!             nconpt(ihpa),iconpt(1,ihpa),conpt(1,ihpa))
 
       else
         write(6,"('problem 1')")
@@ -1280,7 +1287,7 @@
 
   enddo ! --- by ish
 
-! --- evaluate perturbations in per cent
+! --- evaluate perturbations in percent
 
   dvsh=vsh3drel
   dvsv=vsv3drel
@@ -1331,14 +1338,19 @@
       numcof=numcoe(ihpa)
 
 ! originally called
-!        call splcon(y,x,numcof,xlaspl(1,ihpa), &
-!              xlospl(1,ihpa),radspl(1,ihpa), &
-!              nconpt(ihpa),iconpt(1,ihpa),conpt(1,ihpa))
+         call splcon(y,x,numcof,xlaspl(1,ihpa), &
+               xlospl(1,ihpa),radspl(1,ihpa), &
+               nconpt(ihpa),iconpt(1,ihpa),conpt(1,ihpa))
 
-! making sure array bounds
-      call splcon(y,x,numcof,xlaspl(1:numcof,ihpa), &
-              xlospl(1:numcof,ihpa),radspl(1:numcof,ihpa), &
-              nconpt(ihpa),iconpt(1:maxver,ihpa),conpt(1:maxver,ihpa))
+! making sure of array bounds
+!! note from DK DK, July 2013: the code below is correct but it is exactly the same as the code
+!! note from DK DK, July 2013: above because xlaspl(1:numcof,ihpa) is contiguous in memory, and thus
+!! note from DK DK, July 2013: using a pointer to it such as xlaspl(1,ihpa) in the call seems fine as well;
+!! note from DK DK, July 2013: and some compilers could create expensive and useless memory copies for each call
+!! note from DK DK, July 2013: with the new syntax below
+!       call splcon(y,x,numcof,xlaspl(1:numcof,ihpa), &
+!             xlospl(1:numcof,ihpa),radspl(1:numcof,ihpa), &
+!             nconpt(ihpa),iconpt(1:maxver,ihpa),conpt(1:maxver,ihpa))
 
 
     else

Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/model_s40rts.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/model_s40rts.f90	2013-07-24 14:30:04 UTC (rev 22661)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/model_s40rts.f90	2013-07-24 15:33:12 UTC (rev 22662)
@@ -67,13 +67,12 @@
 
 ! standard routine to setup model
 
+  use mpi
   use model_s40rts_par
 
   implicit none
 
   include "constants.h"
-  ! standard include of the MPI library
-  include 'mpif.h'
 
   integer :: myrank
 

Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/rules.mk
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/rules.mk	2013-07-24 14:30:04 UTC (rev 22661)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/rules.mk	2013-07-24 15:33:12 UTC (rev 22662)
@@ -104,7 +104,7 @@
 # These files come from the shared directory
 meshfem3D_SHARED_OBJECTS = \
 	$O/auto_ner.o \
-	$O/broadcast_compute_parameters.o \
+	$O/broadcast_computed_parameters.o \
 	$O/calendar.o \
 	$O/count_number_of_sources.o \
 	$O/create_name_database.o \

Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/save_arrays_solver.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/save_arrays_solver.f90	2013-07-24 14:30:04 UTC (rev 22661)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/save_arrays_solver.f90	2013-07-24 15:33:12 UTC (rev 22662)
@@ -237,7 +237,7 @@
 
   ! mass matrices
   !
-  ! in the case of stacey boundary conditions, add C*deltat/2 contribution to the mass matrix
+  ! in the case of Stacey boundary conditions, add C*deltat/2 contribution to the mass matrix
   ! on the Stacey edges for the crust_mantle and outer_core regions but not for the inner_core region
   ! thus the mass matrix must be replaced by three mass matrices including the "C" damping matrix
   !
@@ -294,7 +294,6 @@
     open(unit=27, file=prname(1:len_trim(prname))//'attenuation.bin', &
           status='unknown', form='unformatted',action='write',iostat=ier)
     if( ier /= 0 ) call exit_mpi(myrank,'error opening attenuation.bin file')
-
     write(27) tau_s
     write(27) tau_e_store
     write(27) Qmu_store

Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/model_1066a.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/model_1066a.f90	2013-07-24 14:30:04 UTC (rev 22661)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/model_1066a.f90	2013-07-24 15:33:12 UTC (rev 22662)
@@ -144,9 +144,9 @@
 ! make sure Vs is zero in the outer core even if roundoff errors on depth
 ! also set fictitious attenuation to a very high value (attenuation is not used in the fluid)
   if(iregion_code == IREGION_OUTER_CORE) then
-    vs = 0.d0
-    Qkappa = 3000.d0
-    Qmu = 3000.d0
+    vs = 0.dd0
+    Qkappa = 3000.dd0
+    Qmu = 3000.dd0
   endif
 
 ! non-dimensionalize
@@ -181,653 +181,653 @@
 
   logical USE_EXTERNAL_CRUSTAL_MODEL
 
-  integer i
+  integer :: i
 
 ! define all the values in the model
 
-  M1066a_V%radius_1066a(  1) =  0.000000000000000
-  M1066a_V%radius_1066a(  2) =   38400.0000000000
-  M1066a_V%radius_1066a(  3) =   76810.0000000000
-  M1066a_V%radius_1066a(  4) =   115210.000000000
-  M1066a_V%radius_1066a(  5) =   153610.000000000
-  M1066a_V%radius_1066a(  6) =   192020.000000000
-  M1066a_V%radius_1066a(  7) =   230420.000000000
-  M1066a_V%radius_1066a(  8) =   268820.000000000
-  M1066a_V%radius_1066a(  9) =   307220.000000000
-  M1066a_V%radius_1066a( 10) =   345630.000000000
-  M1066a_V%radius_1066a( 11) =   384030.000000000
-  M1066a_V%radius_1066a( 12) =   422430.000000000
-  M1066a_V%radius_1066a( 13) =   460840.000000000
-  M1066a_V%radius_1066a( 14) =   499240.000000000
-  M1066a_V%radius_1066a( 15) =   537640.000000000
-  M1066a_V%radius_1066a( 16) =   576050.000000000
-  M1066a_V%radius_1066a( 17) =   614450.000000000
-  M1066a_V%radius_1066a( 18) =   652850.000000000
-  M1066a_V%radius_1066a( 19) =   691260.000000000
-  M1066a_V%radius_1066a( 20) =   729660.000000000
-  M1066a_V%radius_1066a( 21) =   768060.000000000
-  M1066a_V%radius_1066a( 22) =   806460.000000000
-  M1066a_V%radius_1066a( 23) =   844870.000000000
-  M1066a_V%radius_1066a( 24) =   883270.000000000
-  M1066a_V%radius_1066a( 25) =   921670.000000000
-  M1066a_V%radius_1066a( 26) =   960080.000000000
-  M1066a_V%radius_1066a( 27) =   998480.000000000
-  M1066a_V%radius_1066a( 28) =   1036880.00000000
-  M1066a_V%radius_1066a( 29) =   1075290.00000000
-  M1066a_V%radius_1066a( 30) =   1113690.00000000
-  M1066a_V%radius_1066a( 31) =   1152090.00000000
-  M1066a_V%radius_1066a( 32) =   1190500.00000000
-  M1066a_V%radius_1066a( 33) =   1229480.00000000
-  M1066a_V%radius_1066a( 34) =   1229480.00000000
-  M1066a_V%radius_1066a( 35) =   1299360.00000000
-  M1066a_V%radius_1066a( 36) =   1369820.00000000
-  M1066a_V%radius_1066a( 37) =   1440280.00000000
-  M1066a_V%radius_1066a( 38) =   1510740.00000000
-  M1066a_V%radius_1066a( 39) =   1581190.00000000
-  M1066a_V%radius_1066a( 40) =   1651650.00000000
-  M1066a_V%radius_1066a( 41) =   1722110.00000000
-  M1066a_V%radius_1066a( 42) =   1792570.00000000
-  M1066a_V%radius_1066a( 43) =   1863030.00000000
-  M1066a_V%radius_1066a( 44) =   1933490.00000000
-  M1066a_V%radius_1066a( 45) =   2003950.00000000
-  M1066a_V%radius_1066a( 46) =   2074410.00000000
-  M1066a_V%radius_1066a( 47) =   2144870.00000000
-  M1066a_V%radius_1066a( 48) =   2215330.00000000
-  M1066a_V%radius_1066a( 49) =   2285790.00000000
-  M1066a_V%radius_1066a( 50) =   2356240.00000000
-  M1066a_V%radius_1066a( 51) =   2426700.00000000
-  M1066a_V%radius_1066a( 52) =   2497160.00000000
-  M1066a_V%radius_1066a( 53) =   2567620.00000000
-  M1066a_V%radius_1066a( 54) =   2638080.00000000
-  M1066a_V%radius_1066a( 55) =   2708540.00000000
-  M1066a_V%radius_1066a( 56) =   2779000.00000000
-  M1066a_V%radius_1066a( 57) =   2849460.00000000
-  M1066a_V%radius_1066a( 58) =   2919920.00000000
-  M1066a_V%radius_1066a( 59) =   2990380.00000000
-  M1066a_V%radius_1066a( 60) =   3060840.00000000
-  M1066a_V%radius_1066a( 61) =   3131300.00000000
-  M1066a_V%radius_1066a( 62) =   3201750.00000000
-  M1066a_V%radius_1066a( 63) =   3272210.00000000
-  M1066a_V%radius_1066a( 64) =   3342670.00000000
-  M1066a_V%radius_1066a( 65) =   3413130.00000000
-  M1066a_V%radius_1066a( 66) =   3484300.00000000
-  M1066a_V%radius_1066a( 67) =   3484300.00000000
-  M1066a_V%radius_1066a( 68) =   3518220.00000000
-  M1066a_V%radius_1066a( 69) =   3552850.00000000
-  M1066a_V%radius_1066a( 70) =   3587490.00000000
-  M1066a_V%radius_1066a( 71) =   3622120.00000000
-  M1066a_V%radius_1066a( 72) =   3656750.00000000
-  M1066a_V%radius_1066a( 73) =   3691380.00000000
-  M1066a_V%radius_1066a( 74) =   3726010.00000000
-  M1066a_V%radius_1066a( 75) =   3760640.00000000
-  M1066a_V%radius_1066a( 76) =   3795270.00000000
-  M1066a_V%radius_1066a( 77) =   3829910.00000000
-  M1066a_V%radius_1066a( 78) =   3864540.00000000
-  M1066a_V%radius_1066a( 79) =   3899170.00000000
-  M1066a_V%radius_1066a( 80) =   3933800.00000000
-  M1066a_V%radius_1066a( 81) =   3968430.00000000
-  M1066a_V%radius_1066a( 82) =   4003060.00000000
-  M1066a_V%radius_1066a( 83) =   4037690.00000000
-  M1066a_V%radius_1066a( 84) =   4072330.00000000
-  M1066a_V%radius_1066a( 85) =   4106960.00000000
-  M1066a_V%radius_1066a( 86) =   4141590.00000000
-  M1066a_V%radius_1066a( 87) =   4176220.00000000
-  M1066a_V%radius_1066a( 88) =   4210850.00000000
-  M1066a_V%radius_1066a( 89) =   4245480.00000000
-  M1066a_V%radius_1066a( 90) =   4280110.00000000
-  M1066a_V%radius_1066a( 91) =   4314740.00000000
-  M1066a_V%radius_1066a( 92) =   4349380.00000000
-  M1066a_V%radius_1066a( 93) =   4384010.00000000
-  M1066a_V%radius_1066a( 94) =   4418640.00000000
-  M1066a_V%radius_1066a( 95) =   4453270.00000000
-  M1066a_V%radius_1066a( 96) =   4487900.00000000
-  M1066a_V%radius_1066a( 97) =   4522530.00000000
-  M1066a_V%radius_1066a( 98) =   4557160.00000000
-  M1066a_V%radius_1066a( 99) =   4591800.00000000
-  M1066a_V%radius_1066a(100) =   4626430.00000000
-  M1066a_V%radius_1066a(101) =   4661060.00000000
-  M1066a_V%radius_1066a(102) =   4695690.00000000
-  M1066a_V%radius_1066a(103) =   4730320.00000000
-  M1066a_V%radius_1066a(104) =   4764950.00000000
-  M1066a_V%radius_1066a(105) =   4799580.00000000
-  M1066a_V%radius_1066a(106) =   4834220.00000000
-  M1066a_V%radius_1066a(107) =   4868850.00000000
-  M1066a_V%radius_1066a(108) =   4903480.00000000
-  M1066a_V%radius_1066a(109) =   4938110.00000000
-  M1066a_V%radius_1066a(110) =   4972740.00000000
-  M1066a_V%radius_1066a(111) =   5007370.00000000
-  M1066a_V%radius_1066a(112) =   5042000.00000000
-  M1066a_V%radius_1066a(113) =   5076640.00000000
-  M1066a_V%radius_1066a(114) =   5111270.00000000
-  M1066a_V%radius_1066a(115) =   5145900.00000000
-  M1066a_V%radius_1066a(116) =   5180530.00000000
-  M1066a_V%radius_1066a(117) =   5215160.00000000
-  M1066a_V%radius_1066a(118) =   5249790.00000000
-  M1066a_V%radius_1066a(119) =   5284420.00000000
-  M1066a_V%radius_1066a(120) =   5319060.00000000
-  M1066a_V%radius_1066a(121) =   5353690.00000000
-  M1066a_V%radius_1066a(122) =   5388320.00000000
-  M1066a_V%radius_1066a(123) =   5422950.00000000
-  M1066a_V%radius_1066a(124) =   5457580.00000000
-  M1066a_V%radius_1066a(125) =   5492210.00000000
-  M1066a_V%radius_1066a(126) =   5526840.00000000
-  M1066a_V%radius_1066a(127) =   5561470.00000000
-  M1066a_V%radius_1066a(128) =   5596110.00000000
-  M1066a_V%radius_1066a(129) =   5630740.00000000
-  M1066a_V%radius_1066a(130) =   5665370.00000000
-  M1066a_V%radius_1066a(131) =   5700000.00000000
-  M1066a_V%radius_1066a(132) =   5700000.00000000
-  M1066a_V%radius_1066a(133) =   5731250.00000000
-  M1066a_V%radius_1066a(134) =   5762500.00000000
-  M1066a_V%radius_1066a(135) =   5793750.00000000
-  M1066a_V%radius_1066a(136) =   5825000.00000000
-  M1066a_V%radius_1066a(137) =   5856250.00000000
-  M1066a_V%radius_1066a(138) =   5887500.00000000
-  M1066a_V%radius_1066a(139) =   5918750.00000000
-  M1066a_V%radius_1066a(140) =   5950000.00000000
-  M1066a_V%radius_1066a(141) =   5950000.00000000
-  M1066a_V%radius_1066a(142) =   5975630.00000000
-  M1066a_V%radius_1066a(143) =   6001250.00000000
-  M1066a_V%radius_1066a(144) =   6026880.00000000
-  M1066a_V%radius_1066a(145) =   6052500.00000000
-  M1066a_V%radius_1066a(146) =   6078130.00000000
-  M1066a_V%radius_1066a(147) =   6103750.00000000
-  M1066a_V%radius_1066a(148) =   6129380.00000000
-  M1066a_V%radius_1066a(149) =   6155000.00000000
-  M1066a_V%radius_1066a(150) =   6180630.00000000
-  M1066a_V%radius_1066a(151) =   6206250.00000000
-  M1066a_V%radius_1066a(152) =   6231880.00000000
-  M1066a_V%radius_1066a(153) =   6257500.00000000
-  M1066a_V%radius_1066a(154) =   6283130.00000000
-  M1066a_V%radius_1066a(155) =   6308750.00000000
-  M1066a_V%radius_1066a(156) =   6334380.00000000
-  M1066a_V%radius_1066a(157) =   6360000.00000000
-  M1066a_V%radius_1066a(158) =   6360000.00000000
-  M1066a_V%radius_1066a(159) =   6365500.00000000
-  M1066a_V%radius_1066a(160) =   6371000.00000000
+  M1066a_V%radius_1066a(  1) =  0.d0
+  M1066a_V%radius_1066a(  2) =   38400.d0
+  M1066a_V%radius_1066a(  3) =   76810.d0
+  M1066a_V%radius_1066a(  4) =   115210.d0
+  M1066a_V%radius_1066a(  5) =   153610.d0
+  M1066a_V%radius_1066a(  6) =   192020.d0
+  M1066a_V%radius_1066a(  7) =   230420.d0
+  M1066a_V%radius_1066a(  8) =   268820.d0
+  M1066a_V%radius_1066a(  9) =   307220.d0
+  M1066a_V%radius_1066a( 10) =   345630.d0
+  M1066a_V%radius_1066a( 11) =   384030.d0
+  M1066a_V%radius_1066a( 12) =   422430.d0
+  M1066a_V%radius_1066a( 13) =   460840.d0
+  M1066a_V%radius_1066a( 14) =   499240.d0
+  M1066a_V%radius_1066a( 15) =   537640.d0
+  M1066a_V%radius_1066a( 16) =   576050.d0
+  M1066a_V%radius_1066a( 17) =   614450.d0
+  M1066a_V%radius_1066a( 18) =   652850.d0
+  M1066a_V%radius_1066a( 19) =   691260.d0
+  M1066a_V%radius_1066a( 20) =   729660.d0
+  M1066a_V%radius_1066a( 21) =   768060.d0
+  M1066a_V%radius_1066a( 22) =   806460.d0
+  M1066a_V%radius_1066a( 23) =   844870.d0
+  M1066a_V%radius_1066a( 24) =   883270.d0
+  M1066a_V%radius_1066a( 25) =   921670.d0
+  M1066a_V%radius_1066a( 26) =   960080.d0
+  M1066a_V%radius_1066a( 27) =   998480.d0
+  M1066a_V%radius_1066a( 28) =   1036880.d0
+  M1066a_V%radius_1066a( 29) =   1075290.d0
+  M1066a_V%radius_1066a( 30) =   1113690.d0
+  M1066a_V%radius_1066a( 31) =   1152090.d0
+  M1066a_V%radius_1066a( 32) =   1190500.d0
+  M1066a_V%radius_1066a( 33) =   1229480.d0
+  M1066a_V%radius_1066a( 34) =   1229480.d0
+  M1066a_V%radius_1066a( 35) =   1299360.d0
+  M1066a_V%radius_1066a( 36) =   1369820.d0
+  M1066a_V%radius_1066a( 37) =   1440280.d0
+  M1066a_V%radius_1066a( 38) =   1510740.d0
+  M1066a_V%radius_1066a( 39) =   1581190.d0
+  M1066a_V%radius_1066a( 40) =   1651650.d0
+  M1066a_V%radius_1066a( 41) =   1722110.d0
+  M1066a_V%radius_1066a( 42) =   1792570.d0
+  M1066a_V%radius_1066a( 43) =   1863030.d0
+  M1066a_V%radius_1066a( 44) =   1933490.d0
+  M1066a_V%radius_1066a( 45) =   2003950.d0
+  M1066a_V%radius_1066a( 46) =   2074410.d0
+  M1066a_V%radius_1066a( 47) =   2144870.d0
+  M1066a_V%radius_1066a( 48) =   2215330.d0
+  M1066a_V%radius_1066a( 49) =   2285790.d0
+  M1066a_V%radius_1066a( 50) =   2356240.d0
+  M1066a_V%radius_1066a( 51) =   2426700.d0
+  M1066a_V%radius_1066a( 52) =   2497160.d0
+  M1066a_V%radius_1066a( 53) =   2567620.d0
+  M1066a_V%radius_1066a( 54) =   2638080.d0
+  M1066a_V%radius_1066a( 55) =   2708540.d0
+  M1066a_V%radius_1066a( 56) =   2779000.d0
+  M1066a_V%radius_1066a( 57) =   2849460.d0
+  M1066a_V%radius_1066a( 58) =   2919920.d0
+  M1066a_V%radius_1066a( 59) =   2990380.d0
+  M1066a_V%radius_1066a( 60) =   3060840.d0
+  M1066a_V%radius_1066a( 61) =   3131300.d0
+  M1066a_V%radius_1066a( 62) =   3201750.d0
+  M1066a_V%radius_1066a( 63) =   3272210.d0
+  M1066a_V%radius_1066a( 64) =   3342670.d0
+  M1066a_V%radius_1066a( 65) =   3413130.d0
+  M1066a_V%radius_1066a( 66) =   3484300.d0
+  M1066a_V%radius_1066a( 67) =   3484300.d0
+  M1066a_V%radius_1066a( 68) =   3518220.d0
+  M1066a_V%radius_1066a( 69) =   3552850.d0
+  M1066a_V%radius_1066a( 70) =   3587490.d0
+  M1066a_V%radius_1066a( 71) =   3622120.d0
+  M1066a_V%radius_1066a( 72) =   3656750.d0
+  M1066a_V%radius_1066a( 73) =   3691380.d0
+  M1066a_V%radius_1066a( 74) =   3726010.d0
+  M1066a_V%radius_1066a( 75) =   3760640.d0
+  M1066a_V%radius_1066a( 76) =   3795270.d0
+  M1066a_V%radius_1066a( 77) =   3829910.d0
+  M1066a_V%radius_1066a( 78) =   3864540.d0
+  M1066a_V%radius_1066a( 79) =   3899170.d0
+  M1066a_V%radius_1066a( 80) =   3933800.d0
+  M1066a_V%radius_1066a( 81) =   3968430.d0
+  M1066a_V%radius_1066a( 82) =   4003060.d0
+  M1066a_V%radius_1066a( 83) =   4037690.d0
+  M1066a_V%radius_1066a( 84) =   4072330.d0
+  M1066a_V%radius_1066a( 85) =   4106960.d0
+  M1066a_V%radius_1066a( 86) =   4141590.d0
+  M1066a_V%radius_1066a( 87) =   4176220.d0
+  M1066a_V%radius_1066a( 88) =   4210850.d0
+  M1066a_V%radius_1066a( 89) =   4245480.d0
+  M1066a_V%radius_1066a( 90) =   4280110.d0
+  M1066a_V%radius_1066a( 91) =   4314740.d0
+  M1066a_V%radius_1066a( 92) =   4349380.d0
+  M1066a_V%radius_1066a( 93) =   4384010.d0
+  M1066a_V%radius_1066a( 94) =   4418640.d0
+  M1066a_V%radius_1066a( 95) =   4453270.d0
+  M1066a_V%radius_1066a( 96) =   4487900.d0
+  M1066a_V%radius_1066a( 97) =   4522530.d0
+  M1066a_V%radius_1066a( 98) =   4557160.d0
+  M1066a_V%radius_1066a( 99) =   4591800.d0
+  M1066a_V%radius_1066a(100) =   4626430.d0
+  M1066a_V%radius_1066a(101) =   4661060.d0
+  M1066a_V%radius_1066a(102) =   4695690.d0
+  M1066a_V%radius_1066a(103) =   4730320.d0
+  M1066a_V%radius_1066a(104) =   4764950.d0
+  M1066a_V%radius_1066a(105) =   4799580.d0
+  M1066a_V%radius_1066a(106) =   4834220.d0
+  M1066a_V%radius_1066a(107) =   4868850.d0
+  M1066a_V%radius_1066a(108) =   4903480.d0
+  M1066a_V%radius_1066a(109) =   4938110.d0
+  M1066a_V%radius_1066a(110) =   4972740.d0
+  M1066a_V%radius_1066a(111) =   5007370.d0
+  M1066a_V%radius_1066a(112) =   5042000.d0
+  M1066a_V%radius_1066a(113) =   5076640.d0
+  M1066a_V%radius_1066a(114) =   5111270.d0
+  M1066a_V%radius_1066a(115) =   5145900.d0
+  M1066a_V%radius_1066a(116) =   5180530.d0
+  M1066a_V%radius_1066a(117) =   5215160.d0
+  M1066a_V%radius_1066a(118) =   5249790.d0
+  M1066a_V%radius_1066a(119) =   5284420.d0
+  M1066a_V%radius_1066a(120) =   5319060.d0
+  M1066a_V%radius_1066a(121) =   5353690.d0
+  M1066a_V%radius_1066a(122) =   5388320.d0
+  M1066a_V%radius_1066a(123) =   5422950.d0
+  M1066a_V%radius_1066a(124) =   5457580.d0
+  M1066a_V%radius_1066a(125) =   5492210.d0
+  M1066a_V%radius_1066a(126) =   5526840.d0
+  M1066a_V%radius_1066a(127) =   5561470.d0
+  M1066a_V%radius_1066a(128) =   5596110.d0
+  M1066a_V%radius_1066a(129) =   5630740.d0
+  M1066a_V%radius_1066a(130) =   5665370.d0
+  M1066a_V%radius_1066a(131) =   5700000.d0
+  M1066a_V%radius_1066a(132) =   5700000.d0
+  M1066a_V%radius_1066a(133) =   5731250.d0
+  M1066a_V%radius_1066a(134) =   5762500.d0
+  M1066a_V%radius_1066a(135) =   5793750.d0
+  M1066a_V%radius_1066a(136) =   5825000.d0
+  M1066a_V%radius_1066a(137) =   5856250.d0
+  M1066a_V%radius_1066a(138) =   5887500.d0
+  M1066a_V%radius_1066a(139) =   5918750.d0
+  M1066a_V%radius_1066a(140) =   5950000.d0
+  M1066a_V%radius_1066a(141) =   5950000.d0
+  M1066a_V%radius_1066a(142) =   5975630.d0
+  M1066a_V%radius_1066a(143) =   6001250.d0
+  M1066a_V%radius_1066a(144) =   6026880.d0
+  M1066a_V%radius_1066a(145) =   6052500.d0
+  M1066a_V%radius_1066a(146) =   6078130.d0
+  M1066a_V%radius_1066a(147) =   6103750.d0
+  M1066a_V%radius_1066a(148) =   6129380.d0
+  M1066a_V%radius_1066a(149) =   6155000.d0
+  M1066a_V%radius_1066a(150) =   6180630.d0
+  M1066a_V%radius_1066a(151) =   6206250.d0
+  M1066a_V%radius_1066a(152) =   6231880.d0
+  M1066a_V%radius_1066a(153) =   6257500.d0
+  M1066a_V%radius_1066a(154) =   6283130.d0
+  M1066a_V%radius_1066a(155) =   6308750.d0
+  M1066a_V%radius_1066a(156) =   6334380.d0
+  M1066a_V%radius_1066a(157) =   6360000.d0
+  M1066a_V%radius_1066a(158) =   6360000.d0
+  M1066a_V%radius_1066a(159) =   6365500.d0
+  M1066a_V%radius_1066a(160) =   6371000.d0
 
-  M1066a_V%density_1066a(  1) =   13.4290300000000
-  M1066a_V%density_1066a(  2) =   13.4256300000000
-  M1066a_V%density_1066a(  3) =   13.4191300000000
-  M1066a_V%density_1066a(  4) =   13.4135300000000
-  M1066a_V%density_1066a(  5) =   13.4072300000000
-  M1066a_V%density_1066a(  6) =   13.4003200000000
-  M1066a_V%density_1066a(  7) =   13.3929200000000
-  M1066a_V%density_1066a(  8) =   13.3847100000000
-  M1066a_V%density_1066a(  9) =   13.3754000000000
-  M1066a_V%density_1066a( 10) =   13.3649000000000
-  M1066a_V%density_1066a( 11) =   13.3527900000000
-  M1066a_V%density_1066a( 12) =   13.3389800000000
-  M1066a_V%density_1066a( 13) =   13.3238700000000
-  M1066a_V%density_1066a( 14) =   13.3078500000000
-  M1066a_V%density_1066a( 15) =   13.2914400000000
-  M1066a_V%density_1066a( 16) =   13.2750300000000
-  M1066a_V%density_1066a( 17) =   13.2589100000000
-  M1066a_V%density_1066a( 18) =   13.2431000000000
-  M1066a_V%density_1066a( 19) =   13.2275800000000
-  M1066a_V%density_1066a( 20) =   13.2123600000000
-  M1066a_V%density_1066a( 21) =   13.1972500000000
-  M1066a_V%density_1066a( 22) =   13.1823300000000
-  M1066a_V%density_1066a( 23) =   13.1675100000000
-  M1066a_V%density_1066a( 24) =   13.1527800000000
-  M1066a_V%density_1066a( 25) =   13.1382600000000
-  M1066a_V%density_1066a( 26) =   13.1239400000000
-  M1066a_V%density_1066a( 27) =   13.1095200000000
-  M1066a_V%density_1066a( 28) =   13.0953900000000
-  M1066a_V%density_1066a( 29) =   13.0811600000000
-  M1066a_V%density_1066a( 30) =   13.0670400000000
-  M1066a_V%density_1066a( 31) =   13.0525100000000
-  M1066a_V%density_1066a( 32) =   13.0385800000000
-  M1066a_V%density_1066a( 33) =   13.0287500000000
-  M1066a_V%density_1066a( 34) =   12.1606500000000
-  M1066a_V%density_1066a( 35) =   12.1169900000000
-  M1066a_V%density_1066a( 36) =   12.0748300000000
-  M1066a_V%density_1066a( 37) =   12.0330700000000
-  M1066a_V%density_1066a( 38) =   11.9916000000000
-  M1066a_V%density_1066a( 39) =   11.9507300000000
-  M1066a_V%density_1066a( 40) =   11.9104600000000
-  M1066a_V%density_1066a( 41) =   11.8693800000000
-  M1066a_V%density_1066a( 42) =   11.8248100000000
-  M1066a_V%density_1066a( 43) =   11.7753200000000
-  M1066a_V%density_1066a( 44) =   11.7220400000000
-  M1066a_V%density_1066a( 45) =   11.6665500000000
-  M1066a_V%density_1066a( 46) =   11.6085600000000
-  M1066a_V%density_1066a( 47) =   11.5469600000000
-  M1066a_V%density_1066a( 48) =   11.4809600000000
-  M1066a_V%density_1066a( 49) =   11.4116600000000
-  M1066a_V%density_1066a( 50) =   11.3411600000000
-  M1066a_V%density_1066a( 51) =   11.2705500000000
-  M1066a_V%density_1066a( 52) =   11.1982400000000
-  M1066a_V%density_1066a( 53) =   11.1214200000000
-  M1066a_V%density_1066a( 54) =   11.0384100000000
-  M1066a_V%density_1066a( 55) =   10.9511900000000
-  M1066a_V%density_1066a( 56) =   10.8631600000000
-  M1066a_V%density_1066a( 57) =   10.7770300000000
-  M1066a_V%density_1066a( 58) =   10.6925000000000
-  M1066a_V%density_1066a( 59) =   10.6076700000000
-  M1066a_V%density_1066a( 60) =   10.5207300000000
-  M1066a_V%density_1066a( 61) =   10.4312000000000
-  M1066a_V%density_1066a( 62) =   10.3377500000000
-  M1066a_V%density_1066a( 63) =   10.2396100000000
-  M1066a_V%density_1066a( 64) =   10.1378600000000
-  M1066a_V%density_1066a( 65) =   10.0323000000000
-  M1066a_V%density_1066a( 66) =   9.91745000000000
-  M1066a_V%density_1066a( 67) =   5.53205000000000
-  M1066a_V%density_1066a( 68) =   5.52147000000000
-  M1066a_V%density_1066a( 69) =   5.50959000000000
-  M1066a_V%density_1066a( 70) =   5.49821000000000
-  M1066a_V%density_1066a( 71) =   5.48673000000000
-  M1066a_V%density_1066a( 72) =   5.47495000000000
-  M1066a_V%density_1066a( 73) =   5.46297000000000
-  M1066a_V%density_1066a( 74) =   5.45049000000000
-  M1066a_V%density_1066a( 75) =   5.43741000000000
-  M1066a_V%density_1066a( 76) =   5.42382000000000
-  M1066a_V%density_1066a( 77) =   5.40934000000000
-  M1066a_V%density_1066a( 78) =   5.39375000000000
-  M1066a_V%density_1066a( 79) =   5.37717000000000
-  M1066a_V%density_1066a( 80) =   5.35958000000000
-  M1066a_V%density_1066a( 81) =   5.34079000000000
-  M1066a_V%density_1066a( 82) =   5.32100000000000
-  M1066a_V%density_1066a( 83) =   5.30031000000000
-  M1066a_V%density_1066a( 84) =   5.27902000000000
-  M1066a_V%density_1066a( 85) =   5.25733000000000
-  M1066a_V%density_1066a( 86) =   5.23554000000000
-  M1066a_V%density_1066a( 87) =   5.21375000000000
-  M1066a_V%density_1066a( 88) =   5.19196000000000
-  M1066a_V%density_1066a( 89) =   5.17056000000000
-  M1066a_V%density_1066a( 90) =   5.14937000000000
-  M1066a_V%density_1066a( 91) =   5.12827000000000
-  M1066a_V%density_1066a( 92) =   5.10758000000000
-  M1066a_V%density_1066a( 93) =   5.08728000000000
-  M1066a_V%density_1066a( 94) =   5.06738000000000
-  M1066a_V%density_1066a( 95) =   5.04769000000000
-  M1066a_V%density_1066a( 96) =   5.02809000000000
-  M1066a_V%density_1066a( 97) =   5.00869000000000
-  M1066a_V%density_1066a( 98) =   4.98929000000000
-  M1066a_V%density_1066a( 99) =   4.96968000000000
-  M1066a_V%density_1066a(100) =   4.95008000000000
-  M1066a_V%density_1066a(101) =   4.93048000000000
-  M1066a_V%density_1066a(102) =   4.91128000000000
-  M1066a_V%density_1066a(103) =   4.89257000000000
-  M1066a_V%density_1066a(104) =   4.87447000000000
-  M1066a_V%density_1066a(105) =   4.85716000000000
-  M1066a_V%density_1066a(106) =   4.84095000000000
-  M1066a_V%density_1066a(107) =   4.82554000000000
-  M1066a_V%density_1066a(108) =   4.81084000000000
-  M1066a_V%density_1066a(109) =   4.79683000000000
-  M1066a_V%density_1066a(110) =   4.78312000000000
-  M1066a_V%density_1066a(111) =   4.76951000000000
-  M1066a_V%density_1066a(112) =   4.75530000000000
-  M1066a_V%density_1066a(113) =   4.74008000000000
-  M1066a_V%density_1066a(114) =   4.72317000000000
-  M1066a_V%density_1066a(115) =   4.70426000000000
-  M1066a_V%density_1066a(116) =   4.68264000000000
-  M1066a_V%density_1066a(117) =   4.65863000000000
-  M1066a_V%density_1066a(118) =   4.63351000000000
-  M1066a_V%density_1066a(119) =   4.60859000000000
-  M1066a_V%density_1066a(120) =   4.58538000000000
-  M1066a_V%density_1066a(121) =   4.56536000000000
-  M1066a_V%density_1066a(122) =   4.55044000000000
-  M1066a_V%density_1066a(123) =   4.54072000000000
-  M1066a_V%density_1066a(124) =   4.53480000000000
-  M1066a_V%density_1066a(125) =   4.53478000000000
-  M1066a_V%density_1066a(126) =   4.53275000000000
-  M1066a_V%density_1066a(127) =   4.50893000000000
-  M1066a_V%density_1066a(128) =   4.46541000000000
-  M1066a_V%density_1066a(129) =   4.40098000000000
-  M1066a_V%density_1066a(130) =   4.31686000000000
-  M1066a_V%density_1066a(131) =   4.20553000000000
-  M1066a_V%density_1066a(132) =   4.20553000000000
-  M1066a_V%density_1066a(133) =   4.10272000000000
-  M1066a_V%density_1066a(134) =   4.02250000000000
-  M1066a_V%density_1066a(135) =   3.95789000000000
-  M1066a_V%density_1066a(136) =   3.89997000000000
-  M1066a_V%density_1066a(137) =   3.84675000000000
-  M1066a_V%density_1066a(138) =   3.80144000000000
-  M1066a_V%density_1066a(139) =   3.76072000000000
-  M1066a_V%density_1066a(140) =   3.70840000000000
-  M1066a_V%density_1066a(141) =   3.70840000000000
-  M1066a_V%density_1066a(142) =   3.65370000000000
-  M1066a_V%density_1066a(143) =   3.59640000000000
-  M1066a_V%density_1066a(144) =   3.54731000000000
-  M1066a_V%density_1066a(145) =   3.50511000000000
-  M1066a_V%density_1066a(146) =   3.46861000000000
-  M1066a_V%density_1066a(147) =   3.43851000000000
-  M1066a_V%density_1066a(148) =   3.41471000000000
-  M1066a_V%density_1066a(149) =   3.39751000000000
-  M1066a_V%density_1066a(150) =   3.38820000000000
-  M1066a_V%density_1066a(151) =   3.38200000000000
-  M1066a_V%density_1066a(152) =   3.37450000000000
-  M1066a_V%density_1066a(153) =   3.36710000000000
-  M1066a_V%density_1066a(154) =   3.35980000000000
-  M1066a_V%density_1066a(155) =   3.35259000000000
-  M1066a_V%density_1066a(156) =   3.34549000000000
-  M1066a_V%density_1066a(157) =   3.33828000000000
-  M1066a_V%density_1066a(158) =   2.17798000000000
-  M1066a_V%density_1066a(159) =   2.17766000000000
-  M1066a_V%density_1066a(160) =   2.17734000000000
+  M1066a_V%density_1066a(  1) =   13.42903d0
+  M1066a_V%density_1066a(  2) =   13.42563d0
+  M1066a_V%density_1066a(  3) =   13.41913d0
+  M1066a_V%density_1066a(  4) =   13.41353d0
+  M1066a_V%density_1066a(  5) =   13.40723d0
+  M1066a_V%density_1066a(  6) =   13.40032d0
+  M1066a_V%density_1066a(  7) =   13.39292d0
+  M1066a_V%density_1066a(  8) =   13.38471d0
+  M1066a_V%density_1066a(  9) =   13.3754d0
+  M1066a_V%density_1066a( 10) =   13.3649d0
+  M1066a_V%density_1066a( 11) =   13.35279d0
+  M1066a_V%density_1066a( 12) =   13.33898d0
+  M1066a_V%density_1066a( 13) =   13.32387d0
+  M1066a_V%density_1066a( 14) =   13.30785d0
+  M1066a_V%density_1066a( 15) =   13.29144d0
+  M1066a_V%density_1066a( 16) =   13.27503d0
+  M1066a_V%density_1066a( 17) =   13.25891d0
+  M1066a_V%density_1066a( 18) =   13.2431d0
+  M1066a_V%density_1066a( 19) =   13.22758d0
+  M1066a_V%density_1066a( 20) =   13.21236d0
+  M1066a_V%density_1066a( 21) =   13.19725d0
+  M1066a_V%density_1066a( 22) =   13.18233d0
+  M1066a_V%density_1066a( 23) =   13.16751d0
+  M1066a_V%density_1066a( 24) =   13.15278d0
+  M1066a_V%density_1066a( 25) =   13.13826d0
+  M1066a_V%density_1066a( 26) =   13.12394d0
+  M1066a_V%density_1066a( 27) =   13.10952d0
+  M1066a_V%density_1066a( 28) =   13.09539d0
+  M1066a_V%density_1066a( 29) =   13.08116d0
+  M1066a_V%density_1066a( 30) =   13.06704d0
+  M1066a_V%density_1066a( 31) =   13.05251d0
+  M1066a_V%density_1066a( 32) =   13.03858d0
+  M1066a_V%density_1066a( 33) =   13.02875d0
+  M1066a_V%density_1066a( 34) =   12.16065d0
+  M1066a_V%density_1066a( 35) =   12.11699d0
+  M1066a_V%density_1066a( 36) =   12.07483d0
+  M1066a_V%density_1066a( 37) =   12.03307d0
+  M1066a_V%density_1066a( 38) =   11.9916d0
+  M1066a_V%density_1066a( 39) =   11.95073d0
+  M1066a_V%density_1066a( 40) =   11.91046d0
+  M1066a_V%density_1066a( 41) =   11.86938d0
+  M1066a_V%density_1066a( 42) =   11.82481d0
+  M1066a_V%density_1066a( 43) =   11.77532d0
+  M1066a_V%density_1066a( 44) =   11.72204d0
+  M1066a_V%density_1066a( 45) =   11.66655d0
+  M1066a_V%density_1066a( 46) =   11.60856d0
+  M1066a_V%density_1066a( 47) =   11.54696d0
+  M1066a_V%density_1066a( 48) =   11.48096d0
+  M1066a_V%density_1066a( 49) =   11.41166d0
+  M1066a_V%density_1066a( 50) =   11.34116d0
+  M1066a_V%density_1066a( 51) =   11.27055d0
+  M1066a_V%density_1066a( 52) =   11.19824d0
+  M1066a_V%density_1066a( 53) =   11.12142d0
+  M1066a_V%density_1066a( 54) =   11.03841d0
+  M1066a_V%density_1066a( 55) =   10.95119d0
+  M1066a_V%density_1066a( 56) =   10.86316d0
+  M1066a_V%density_1066a( 57) =   10.77703d0
+  M1066a_V%density_1066a( 58) =   10.6925d0
+  M1066a_V%density_1066a( 59) =   10.60767d0
+  M1066a_V%density_1066a( 60) =   10.52073d0
+  M1066a_V%density_1066a( 61) =   10.4312d0
+  M1066a_V%density_1066a( 62) =   10.33775d0
+  M1066a_V%density_1066a( 63) =   10.23961d0
+  M1066a_V%density_1066a( 64) =   10.13786d0
+  M1066a_V%density_1066a( 65) =   10.0323d0
+  M1066a_V%density_1066a( 66) =   9.91745d0
+  M1066a_V%density_1066a( 67) =   5.53205d0
+  M1066a_V%density_1066a( 68) =   5.52147d0
+  M1066a_V%density_1066a( 69) =   5.50959d0
+  M1066a_V%density_1066a( 70) =   5.49821d0
+  M1066a_V%density_1066a( 71) =   5.48673d0
+  M1066a_V%density_1066a( 72) =   5.47495d0
+  M1066a_V%density_1066a( 73) =   5.46297d0
+  M1066a_V%density_1066a( 74) =   5.45049d0
+  M1066a_V%density_1066a( 75) =   5.43741d0
+  M1066a_V%density_1066a( 76) =   5.42382d0
+  M1066a_V%density_1066a( 77) =   5.40934d0
+  M1066a_V%density_1066a( 78) =   5.39375d0
+  M1066a_V%density_1066a( 79) =   5.37717d0
+  M1066a_V%density_1066a( 80) =   5.35958d0
+  M1066a_V%density_1066a( 81) =   5.34079d0
+  M1066a_V%density_1066a( 82) =   5.321d0
+  M1066a_V%density_1066a( 83) =   5.30031d0
+  M1066a_V%density_1066a( 84) =   5.27902d0
+  M1066a_V%density_1066a( 85) =   5.25733d0
+  M1066a_V%density_1066a( 86) =   5.23554d0
+  M1066a_V%density_1066a( 87) =   5.21375d0
+  M1066a_V%density_1066a( 88) =   5.19196d0
+  M1066a_V%density_1066a( 89) =   5.17056d0
+  M1066a_V%density_1066a( 90) =   5.14937d0
+  M1066a_V%density_1066a( 91) =   5.12827d0
+  M1066a_V%density_1066a( 92) =   5.10758d0
+  M1066a_V%density_1066a( 93) =   5.08728d0
+  M1066a_V%density_1066a( 94) =   5.06738d0
+  M1066a_V%density_1066a( 95) =   5.04769d0
+  M1066a_V%density_1066a( 96) =   5.02809d0
+  M1066a_V%density_1066a( 97) =   5.00869d0
+  M1066a_V%density_1066a( 98) =   4.98929d0
+  M1066a_V%density_1066a( 99) =   4.96968d0
+  M1066a_V%density_1066a(100) =   4.95008d0
+  M1066a_V%density_1066a(101) =   4.93048d0
+  M1066a_V%density_1066a(102) =   4.91128d0
+  M1066a_V%density_1066a(103) =   4.89257d0
+  M1066a_V%density_1066a(104) =   4.87447d0
+  M1066a_V%density_1066a(105) =   4.85716d0
+  M1066a_V%density_1066a(106) =   4.84095d0
+  M1066a_V%density_1066a(107) =   4.82554d0
+  M1066a_V%density_1066a(108) =   4.81084d0
+  M1066a_V%density_1066a(109) =   4.79683d0
+  M1066a_V%density_1066a(110) =   4.78312d0
+  M1066a_V%density_1066a(111) =   4.76951d0
+  M1066a_V%density_1066a(112) =   4.7553d0
+  M1066a_V%density_1066a(113) =   4.74008d0
+  M1066a_V%density_1066a(114) =   4.72317d0
+  M1066a_V%density_1066a(115) =   4.70426d0
+  M1066a_V%density_1066a(116) =   4.68264d0
+  M1066a_V%density_1066a(117) =   4.65863d0
+  M1066a_V%density_1066a(118) =   4.63351d0
+  M1066a_V%density_1066a(119) =   4.60859d0
+  M1066a_V%density_1066a(120) =   4.58538d0
+  M1066a_V%density_1066a(121) =   4.56536d0
+  M1066a_V%density_1066a(122) =   4.55044d0
+  M1066a_V%density_1066a(123) =   4.54072d0
+  M1066a_V%density_1066a(124) =   4.5348d0
+  M1066a_V%density_1066a(125) =   4.53478d0
+  M1066a_V%density_1066a(126) =   4.53275d0
+  M1066a_V%density_1066a(127) =   4.50893d0
+  M1066a_V%density_1066a(128) =   4.46541d0
+  M1066a_V%density_1066a(129) =   4.40098d0
+  M1066a_V%density_1066a(130) =   4.31686d0
+  M1066a_V%density_1066a(131) =   4.20553d0
+  M1066a_V%density_1066a(132) =   4.20553d0
+  M1066a_V%density_1066a(133) =   4.10272d0
+  M1066a_V%density_1066a(134) =   4.0225d0
+  M1066a_V%density_1066a(135) =   3.95789d0
+  M1066a_V%density_1066a(136) =   3.89997d0
+  M1066a_V%density_1066a(137) =   3.84675d0
+  M1066a_V%density_1066a(138) =   3.80144d0
+  M1066a_V%density_1066a(139) =   3.76072d0
+  M1066a_V%density_1066a(140) =   3.7084d0
+  M1066a_V%density_1066a(141) =   3.7084d0
+  M1066a_V%density_1066a(142) =   3.6537d0
+  M1066a_V%density_1066a(143) =   3.5964d0
+  M1066a_V%density_1066a(144) =   3.54731d0
+  M1066a_V%density_1066a(145) =   3.50511d0
+  M1066a_V%density_1066a(146) =   3.46861d0
+  M1066a_V%density_1066a(147) =   3.43851d0
+  M1066a_V%density_1066a(148) =   3.41471d0
+  M1066a_V%density_1066a(149) =   3.39751d0
+  M1066a_V%density_1066a(150) =   3.3882d0
+  M1066a_V%density_1066a(151) =   3.382d0
+  M1066a_V%density_1066a(152) =   3.3745d0
+  M1066a_V%density_1066a(153) =   3.3671d0
+  M1066a_V%density_1066a(154) =   3.3598d0
+  M1066a_V%density_1066a(155) =   3.35259d0
+  M1066a_V%density_1066a(156) =   3.34549d0
+  M1066a_V%density_1066a(157) =   3.33828d0
+  M1066a_V%density_1066a(158) =   2.17798d0
+  M1066a_V%density_1066a(159) =   2.17766d0
+  M1066a_V%density_1066a(160) =   2.17734d0
 
-  M1066a_V%vp_1066a(  1) =   11.3383000000000
-  M1066a_V%vp_1066a(  2) =   11.3374000000000
-  M1066a_V%vp_1066a(  3) =   11.3347000000000
-  M1066a_V%vp_1066a(  4) =   11.3301000000000
-  M1066a_V%vp_1066a(  5) =   11.3237000000000
-  M1066a_V%vp_1066a(  6) =   11.3155000000000
-  M1066a_V%vp_1066a(  7) =   11.3056000000000
-  M1066a_V%vp_1066a(  8) =   11.2940000000000
-  M1066a_V%vp_1066a(  9) =   11.2810000000000
-  M1066a_V%vp_1066a( 10) =   11.2666000000000
-  M1066a_V%vp_1066a( 11) =   11.2512000000000
-  M1066a_V%vp_1066a( 12) =   11.2349000000000
-  M1066a_V%vp_1066a( 13) =   11.2181000000000
-  M1066a_V%vp_1066a( 14) =   11.2010000000000
-  M1066a_V%vp_1066a( 15) =   11.1840000000000
-  M1066a_V%vp_1066a( 16) =   11.1672000000000
-  M1066a_V%vp_1066a( 17) =   11.1508000000000
-  M1066a_V%vp_1066a( 18) =   11.1351000000000
-  M1066a_V%vp_1066a( 19) =   11.1201000000000
-  M1066a_V%vp_1066a( 20) =   11.1059000000000
-  M1066a_V%vp_1066a( 21) =   11.0924000000000
-  M1066a_V%vp_1066a( 22) =   11.0798000000000
-  M1066a_V%vp_1066a( 23) =   11.0678000000000
-  M1066a_V%vp_1066a( 24) =   11.0564000000000
-  M1066a_V%vp_1066a( 25) =   11.0455000000000
-  M1066a_V%vp_1066a( 26) =   11.0350000000000
-  M1066a_V%vp_1066a( 27) =   11.0248000000000
-  M1066a_V%vp_1066a( 28) =   11.0149000000000
-  M1066a_V%vp_1066a( 29) =   11.0051000000000
-  M1066a_V%vp_1066a( 30) =   10.9953000000000
-  M1066a_V%vp_1066a( 31) =   10.9857000000000
-  M1066a_V%vp_1066a( 32) =   10.9756000000000
-  M1066a_V%vp_1066a( 33) =   10.9687000000000
-  M1066a_V%vp_1066a( 34) =   10.4140000000000
-  M1066a_V%vp_1066a( 35) =   10.3518000000000
-  M1066a_V%vp_1066a( 36) =   10.2922000000000
-  M1066a_V%vp_1066a( 37) =   10.2351000000000
-  M1066a_V%vp_1066a( 38) =   10.1808000000000
-  M1066a_V%vp_1066a( 39) =   10.1297000000000
-  M1066a_V%vp_1066a( 40) =   10.0788000000000
-  M1066a_V%vp_1066a( 41) =   10.0284000000000
-  M1066a_V%vp_1066a( 42) =   9.97880000000000
-  M1066a_V%vp_1066a( 43) =   9.93070000000000
-  M1066a_V%vp_1066a( 44) =   9.88360000000000
-  M1066a_V%vp_1066a( 45) =   9.83530000000000
-  M1066a_V%vp_1066a( 46) =   9.78250000000000
-  M1066a_V%vp_1066a( 47) =   9.72110000000000
-  M1066a_V%vp_1066a( 48) =   9.65210000000000
-  M1066a_V%vp_1066a( 49) =   9.58060000000000
-  M1066a_V%vp_1066a( 50) =   9.51150000000000
-  M1066a_V%vp_1066a( 51) =   9.44650000000000
-  M1066a_V%vp_1066a( 52) =   9.38280000000000
-  M1066a_V%vp_1066a( 53) =   9.31660000000000
-  M1066a_V%vp_1066a( 54) =   9.24420000000000
-  M1066a_V%vp_1066a( 55) =   9.16580000000000
-  M1066a_V%vp_1066a( 56) =   9.08330000000000
-  M1066a_V%vp_1066a( 57) =   8.99870000000000
-  M1066a_V%vp_1066a( 58) =   8.91160000000000
-  M1066a_V%vp_1066a( 59) =   8.82010000000000
-  M1066a_V%vp_1066a( 60) =   8.72230000000000
-  M1066a_V%vp_1066a( 61) =   8.61710000000000
-  M1066a_V%vp_1066a( 62) =   8.50300000000000
-  M1066a_V%vp_1066a( 63) =   8.38070000000000
-  M1066a_V%vp_1066a( 64) =   8.25560000000000
-  M1066a_V%vp_1066a( 65) =   8.13180000000000
-  M1066a_V%vp_1066a( 66) =   8.01120000000000
-  M1066a_V%vp_1066a( 67) =   13.7172000000000
-  M1066a_V%vp_1066a( 68) =   13.7134000000000
-  M1066a_V%vp_1066a( 69) =   13.7089000000000
-  M1066a_V%vp_1066a( 70) =   13.6806000000000
-  M1066a_V%vp_1066a( 71) =   13.6517000000000
-  M1066a_V%vp_1066a( 72) =   13.6251000000000
-  M1066a_V%vp_1066a( 73) =   13.5916000000000
-  M1066a_V%vp_1066a( 74) =   13.5564000000000
-  M1066a_V%vp_1066a( 75) =   13.5165000000000
-  M1066a_V%vp_1066a( 76) =   13.4725000000000
-  M1066a_V%vp_1066a( 77) =   13.4248000000000
-  M1066a_V%vp_1066a( 78) =   13.3742000000000
-  M1066a_V%vp_1066a( 79) =   13.3216000000000
-  M1066a_V%vp_1066a( 80) =   13.2679000000000
-  M1066a_V%vp_1066a( 81) =   13.2142000000000
-  M1066a_V%vp_1066a( 82) =   13.1619000000000
-  M1066a_V%vp_1066a( 83) =   13.1114000000000
-  M1066a_V%vp_1066a( 84) =   13.0631000000000
-  M1066a_V%vp_1066a( 85) =   13.0174000000000
-  M1066a_V%vp_1066a( 86) =   12.9745000000000
-  M1066a_V%vp_1066a( 87) =   12.9346000000000
-  M1066a_V%vp_1066a( 88) =   12.8977000000000
-  M1066a_V%vp_1066a( 89) =   12.8635000000000
-  M1066a_V%vp_1066a( 90) =   12.8318000000000
-  M1066a_V%vp_1066a( 91) =   12.8022000000000
-  M1066a_V%vp_1066a( 92) =   12.7739000000000
-  M1066a_V%vp_1066a( 93) =   12.7463000000000
-  M1066a_V%vp_1066a( 94) =   12.7186000000000
-  M1066a_V%vp_1066a( 95) =   12.6903000000000
-  M1066a_V%vp_1066a( 96) =   12.6610000000000
-  M1066a_V%vp_1066a( 97) =   12.6302000000000
-  M1066a_V%vp_1066a( 98) =   12.5978000000000
-  M1066a_V%vp_1066a( 99) =   12.5637000000000
-  M1066a_V%vp_1066a(100) =   12.5276000000000
-  M1066a_V%vp_1066a(101) =   12.4893000000000
-  M1066a_V%vp_1066a(102) =   12.4485000000000
-  M1066a_V%vp_1066a(103) =   12.4052000000000
-  M1066a_V%vp_1066a(104) =   12.3592000000000
-  M1066a_V%vp_1066a(105) =   12.3105000000000
-  M1066a_V%vp_1066a(106) =   12.2596000000000
-  M1066a_V%vp_1066a(107) =   12.2072000000000
-  M1066a_V%vp_1066a(108) =   12.1538000000000
-  M1066a_V%vp_1066a(109) =   12.0998000000000
-  M1066a_V%vp_1066a(110) =   12.0458000000000
-  M1066a_V%vp_1066a(111) =   11.9920000000000
-  M1066a_V%vp_1066a(112) =   11.9373000000000
-  M1066a_V%vp_1066a(113) =   11.8804000000000
-  M1066a_V%vp_1066a(114) =   11.8200000000000
-  M1066a_V%vp_1066a(115) =   11.7554000000000
-  M1066a_V%vp_1066a(116) =   11.6844000000000
-  M1066a_V%vp_1066a(117) =   11.6079000000000
-  M1066a_V%vp_1066a(118) =   11.5308000000000
-  M1066a_V%vp_1066a(119) =   11.4579000000000
-  M1066a_V%vp_1066a(120) =   11.3935000000000
-  M1066a_V%vp_1066a(121) =   11.3418000000000
-  M1066a_V%vp_1066a(122) =   11.3085000000000
-  M1066a_V%vp_1066a(123) =   11.2938000000000
-  M1066a_V%vp_1066a(124) =   11.2915000000000
-  M1066a_V%vp_1066a(125) =   11.3049000000000
-  M1066a_V%vp_1066a(126) =   11.3123000000000
-  M1066a_V%vp_1066a(127) =   11.2643000000000
-  M1066a_V%vp_1066a(128) =   11.1635000000000
-  M1066a_V%vp_1066a(129) =   11.0063000000000
-  M1066a_V%vp_1066a(130) =   10.7959000000000
-  M1066a_V%vp_1066a(131) =   10.5143000000000
-  M1066a_V%vp_1066a(132) =   10.5143000000000
-  M1066a_V%vp_1066a(133) =   10.2513000000000
-  M1066a_V%vp_1066a(134) =   10.0402000000000
-  M1066a_V%vp_1066a(135) =   9.86480000000000
-  M1066a_V%vp_1066a(136) =   9.70860000000000
-  M1066a_V%vp_1066a(137) =   9.56810000000000
-  M1066a_V%vp_1066a(138) =   9.45120000000000
-  M1066a_V%vp_1066a(139) =   9.35100000000000
-  M1066a_V%vp_1066a(140) =   9.22830000000000
-  M1066a_V%vp_1066a(141) =   9.22830000000000
-  M1066a_V%vp_1066a(142) =   9.10870000000000
-  M1066a_V%vp_1066a(143) =   8.98230000000000
-  M1066a_V%vp_1066a(144) =   8.85920000000000
-  M1066a_V%vp_1066a(145) =   8.73860000000000
-  M1066a_V%vp_1066a(146) =   8.61930000000000
-  M1066a_V%vp_1066a(147) =   8.50180000000000
-  M1066a_V%vp_1066a(148) =   8.38710000000000
-  M1066a_V%vp_1066a(149) =   8.27360000000000
-  M1066a_V%vp_1066a(150) =   8.15850000000000
-  M1066a_V%vp_1066a(151) =   8.05400000000000
-  M1066a_V%vp_1066a(152) =   7.96520000000000
-  M1066a_V%vp_1066a(153) =   7.87340000000000
-  M1066a_V%vp_1066a(154) =   7.79720000000000
-  M1066a_V%vp_1066a(155) =   7.73910000000000
-  M1066a_V%vp_1066a(156) =   7.71340000000000
-  M1066a_V%vp_1066a(157) =   7.70460000000000
-  M1066a_V%vp_1066a(158) =   4.70220000000000
-  M1066a_V%vp_1066a(159) =   4.70010000000000
-  M1066a_V%vp_1066a(160) =   4.69790000000000
+  M1066a_V%vp_1066a(  1) =   11.3383d0
+  M1066a_V%vp_1066a(  2) =   11.3374d0
+  M1066a_V%vp_1066a(  3) =   11.3347d0
+  M1066a_V%vp_1066a(  4) =   11.3301d0
+  M1066a_V%vp_1066a(  5) =   11.3237d0
+  M1066a_V%vp_1066a(  6) =   11.3155d0
+  M1066a_V%vp_1066a(  7) =   11.3056d0
+  M1066a_V%vp_1066a(  8) =   11.294d0
+  M1066a_V%vp_1066a(  9) =   11.281d0
+  M1066a_V%vp_1066a( 10) =   11.2666d0
+  M1066a_V%vp_1066a( 11) =   11.2512d0
+  M1066a_V%vp_1066a( 12) =   11.2349d0
+  M1066a_V%vp_1066a( 13) =   11.2181d0
+  M1066a_V%vp_1066a( 14) =   11.201d0
+  M1066a_V%vp_1066a( 15) =   11.184d0
+  M1066a_V%vp_1066a( 16) =   11.1672d0
+  M1066a_V%vp_1066a( 17) =   11.1508d0
+  M1066a_V%vp_1066a( 18) =   11.1351d0
+  M1066a_V%vp_1066a( 19) =   11.1201d0
+  M1066a_V%vp_1066a( 20) =   11.1059d0
+  M1066a_V%vp_1066a( 21) =   11.0924d0
+  M1066a_V%vp_1066a( 22) =   11.0798d0
+  M1066a_V%vp_1066a( 23) =   11.0678d0
+  M1066a_V%vp_1066a( 24) =   11.0564d0
+  M1066a_V%vp_1066a( 25) =   11.0455d0
+  M1066a_V%vp_1066a( 26) =   11.035d0
+  M1066a_V%vp_1066a( 27) =   11.0248d0
+  M1066a_V%vp_1066a( 28) =   11.0149d0
+  M1066a_V%vp_1066a( 29) =   11.0051d0
+  M1066a_V%vp_1066a( 30) =   10.9953d0
+  M1066a_V%vp_1066a( 31) =   10.9857d0
+  M1066a_V%vp_1066a( 32) =   10.9756d0
+  M1066a_V%vp_1066a( 33) =   10.9687d0
+  M1066a_V%vp_1066a( 34) =   10.414d0
+  M1066a_V%vp_1066a( 35) =   10.3518d0
+  M1066a_V%vp_1066a( 36) =   10.2922d0
+  M1066a_V%vp_1066a( 37) =   10.2351d0
+  M1066a_V%vp_1066a( 38) =   10.1808d0
+  M1066a_V%vp_1066a( 39) =   10.1297d0
+  M1066a_V%vp_1066a( 40) =   10.0788d0
+  M1066a_V%vp_1066a( 41) =   10.0284d0
+  M1066a_V%vp_1066a( 42) =   9.9788d0
+  M1066a_V%vp_1066a( 43) =   9.9307d0
+  M1066a_V%vp_1066a( 44) =   9.8836d0
+  M1066a_V%vp_1066a( 45) =   9.8353d0
+  M1066a_V%vp_1066a( 46) =   9.7825d0
+  M1066a_V%vp_1066a( 47) =   9.7211d0
+  M1066a_V%vp_1066a( 48) =   9.6521d0
+  M1066a_V%vp_1066a( 49) =   9.5806d0
+  M1066a_V%vp_1066a( 50) =   9.5115d0
+  M1066a_V%vp_1066a( 51) =   9.4465d0
+  M1066a_V%vp_1066a( 52) =   9.3828d0
+  M1066a_V%vp_1066a( 53) =   9.3166d0
+  M1066a_V%vp_1066a( 54) =   9.2442d0
+  M1066a_V%vp_1066a( 55) =   9.1658d0
+  M1066a_V%vp_1066a( 56) =   9.0833d0
+  M1066a_V%vp_1066a( 57) =   8.9987d0
+  M1066a_V%vp_1066a( 58) =   8.9116d0
+  M1066a_V%vp_1066a( 59) =   8.8201d0
+  M1066a_V%vp_1066a( 60) =   8.7223d0
+  M1066a_V%vp_1066a( 61) =   8.6171d0
+  M1066a_V%vp_1066a( 62) =   8.503d0
+  M1066a_V%vp_1066a( 63) =   8.3807d0
+  M1066a_V%vp_1066a( 64) =   8.2556d0
+  M1066a_V%vp_1066a( 65) =   8.1318d0
+  M1066a_V%vp_1066a( 66) =   8.0112d0
+  M1066a_V%vp_1066a( 67) =   13.7172d0
+  M1066a_V%vp_1066a( 68) =   13.7134d0
+  M1066a_V%vp_1066a( 69) =   13.7089d0
+  M1066a_V%vp_1066a( 70) =   13.6806d0
+  M1066a_V%vp_1066a( 71) =   13.6517d0
+  M1066a_V%vp_1066a( 72) =   13.6251d0
+  M1066a_V%vp_1066a( 73) =   13.5916d0
+  M1066a_V%vp_1066a( 74) =   13.5564d0
+  M1066a_V%vp_1066a( 75) =   13.5165d0
+  M1066a_V%vp_1066a( 76) =   13.4725d0
+  M1066a_V%vp_1066a( 77) =   13.4248d0
+  M1066a_V%vp_1066a( 78) =   13.3742d0
+  M1066a_V%vp_1066a( 79) =   13.3216d0
+  M1066a_V%vp_1066a( 80) =   13.2679d0
+  M1066a_V%vp_1066a( 81) =   13.2142d0
+  M1066a_V%vp_1066a( 82) =   13.1619d0
+  M1066a_V%vp_1066a( 83) =   13.1114d0
+  M1066a_V%vp_1066a( 84) =   13.0631d0
+  M1066a_V%vp_1066a( 85) =   13.0174d0
+  M1066a_V%vp_1066a( 86) =   12.9745d0
+  M1066a_V%vp_1066a( 87) =   12.9346d0
+  M1066a_V%vp_1066a( 88) =   12.8977d0
+  M1066a_V%vp_1066a( 89) =   12.8635d0
+  M1066a_V%vp_1066a( 90) =   12.8318d0
+  M1066a_V%vp_1066a( 91) =   12.8022d0
+  M1066a_V%vp_1066a( 92) =   12.7739d0
+  M1066a_V%vp_1066a( 93) =   12.7463d0
+  M1066a_V%vp_1066a( 94) =   12.7186d0
+  M1066a_V%vp_1066a( 95) =   12.6903d0
+  M1066a_V%vp_1066a( 96) =   12.661d0
+  M1066a_V%vp_1066a( 97) =   12.6302d0
+  M1066a_V%vp_1066a( 98) =   12.5978d0
+  M1066a_V%vp_1066a( 99) =   12.5637d0
+  M1066a_V%vp_1066a(100) =   12.5276d0
+  M1066a_V%vp_1066a(101) =   12.4893d0
+  M1066a_V%vp_1066a(102) =   12.4485d0
+  M1066a_V%vp_1066a(103) =   12.4052d0
+  M1066a_V%vp_1066a(104) =   12.3592d0
+  M1066a_V%vp_1066a(105) =   12.3105d0
+  M1066a_V%vp_1066a(106) =   12.2596d0
+  M1066a_V%vp_1066a(107) =   12.2072d0
+  M1066a_V%vp_1066a(108) =   12.1538d0
+  M1066a_V%vp_1066a(109) =   12.0998d0
+  M1066a_V%vp_1066a(110) =   12.0458d0
+  M1066a_V%vp_1066a(111) =   11.992d0
+  M1066a_V%vp_1066a(112) =   11.9373d0
+  M1066a_V%vp_1066a(113) =   11.8804d0
+  M1066a_V%vp_1066a(114) =   11.82d0
+  M1066a_V%vp_1066a(115) =   11.7554d0
+  M1066a_V%vp_1066a(116) =   11.6844d0
+  M1066a_V%vp_1066a(117) =   11.6079d0
+  M1066a_V%vp_1066a(118) =   11.5308d0
+  M1066a_V%vp_1066a(119) =   11.4579d0
+  M1066a_V%vp_1066a(120) =   11.3935d0
+  M1066a_V%vp_1066a(121) =   11.3418d0
+  M1066a_V%vp_1066a(122) =   11.3085d0
+  M1066a_V%vp_1066a(123) =   11.2938d0
+  M1066a_V%vp_1066a(124) =   11.2915d0
+  M1066a_V%vp_1066a(125) =   11.3049d0
+  M1066a_V%vp_1066a(126) =   11.3123d0
+  M1066a_V%vp_1066a(127) =   11.2643d0
+  M1066a_V%vp_1066a(128) =   11.1635d0
+  M1066a_V%vp_1066a(129) =   11.0063d0
+  M1066a_V%vp_1066a(130) =   10.7959d0
+  M1066a_V%vp_1066a(131) =   10.5143d0
+  M1066a_V%vp_1066a(132) =   10.5143d0
+  M1066a_V%vp_1066a(133) =   10.2513d0
+  M1066a_V%vp_1066a(134) =   10.0402d0
+  M1066a_V%vp_1066a(135) =   9.8648d0
+  M1066a_V%vp_1066a(136) =   9.7086d0
+  M1066a_V%vp_1066a(137) =   9.5681d0
+  M1066a_V%vp_1066a(138) =   9.4512d0
+  M1066a_V%vp_1066a(139) =   9.351d0
+  M1066a_V%vp_1066a(140) =   9.2283d0
+  M1066a_V%vp_1066a(141) =   9.2283d0
+  M1066a_V%vp_1066a(142) =   9.1087d0
+  M1066a_V%vp_1066a(143) =   8.9823d0
+  M1066a_V%vp_1066a(144) =   8.8592d0
+  M1066a_V%vp_1066a(145) =   8.7386d0
+  M1066a_V%vp_1066a(146) =   8.6193d0
+  M1066a_V%vp_1066a(147) =   8.5018d0
+  M1066a_V%vp_1066a(148) =   8.3871d0
+  M1066a_V%vp_1066a(149) =   8.2736d0
+  M1066a_V%vp_1066a(150) =   8.1585d0
+  M1066a_V%vp_1066a(151) =   8.054d0
+  M1066a_V%vp_1066a(152) =   7.9652d0
+  M1066a_V%vp_1066a(153) =   7.8734d0
+  M1066a_V%vp_1066a(154) =   7.7972d0
+  M1066a_V%vp_1066a(155) =   7.7391d0
+  M1066a_V%vp_1066a(156) =   7.7134d0
+  M1066a_V%vp_1066a(157) =   7.7046d0
+  M1066a_V%vp_1066a(158) =   4.7022d0
+  M1066a_V%vp_1066a(159) =   4.7001d0
+  M1066a_V%vp_1066a(160) =   4.6979d0
 
-  M1066a_V%vs_1066a(  1) =   3.62980000000000
-  M1066a_V%vs_1066a(  2) =   3.62970000000000
-  M1066a_V%vs_1066a(  3) =   3.62940000000000
-  M1066a_V%vs_1066a(  4) =   3.62880000000000
-  M1066a_V%vs_1066a(  5) =   3.62810000000000
-  M1066a_V%vs_1066a(  6) =   3.62710000000000
-  M1066a_V%vs_1066a(  7) =   3.62590000000000
-  M1066a_V%vs_1066a(  8) =   3.62440000000000
-  M1066a_V%vs_1066a(  9) =   3.62280000000000
-  M1066a_V%vs_1066a( 10) =   3.62090000000000
-  M1066a_V%vs_1066a( 11) =   3.61870000000000
-  M1066a_V%vs_1066a( 12) =   3.61630000000000
-  M1066a_V%vs_1066a( 13) =   3.61370000000000
-  M1066a_V%vs_1066a( 14) =   3.61080000000000
-  M1066a_V%vs_1066a( 15) =   3.60760000000000
-  M1066a_V%vs_1066a( 16) =   3.60420000000000
-  M1066a_V%vs_1066a( 17) =   3.60040000000000
-  M1066a_V%vs_1066a( 18) =   3.59650000000000
-  M1066a_V%vs_1066a( 19) =   3.59220000000000
-  M1066a_V%vs_1066a( 20) =   3.58760000000000
-  M1066a_V%vs_1066a( 21) =   3.58280000000000
-  M1066a_V%vs_1066a( 22) =   3.57770000000000
-  M1066a_V%vs_1066a( 23) =   3.57240000000000
-  M1066a_V%vs_1066a( 24) =   3.56680000000000
-  M1066a_V%vs_1066a( 25) =   3.56100000000000
-  M1066a_V%vs_1066a( 26) =   3.55510000000000
-  M1066a_V%vs_1066a( 27) =   3.54900000000000
-  M1066a_V%vs_1066a( 28) =   3.54280000000000
-  M1066a_V%vs_1066a( 29) =   3.53650000000000
-  M1066a_V%vs_1066a( 30) =   3.53010000000000
-  M1066a_V%vs_1066a( 31) =   3.52380000000000
-  M1066a_V%vs_1066a( 32) =   3.51720000000000
-  M1066a_V%vs_1066a( 33) =   3.51180000000000
-  M1066a_V%vs_1066a( 34) =  0.000000000000000
-  M1066a_V%vs_1066a( 35) =  0.000000000000000
-  M1066a_V%vs_1066a( 36) =  0.000000000000000
-  M1066a_V%vs_1066a( 37) =  0.000000000000000
-  M1066a_V%vs_1066a( 38) =  0.000000000000000
-  M1066a_V%vs_1066a( 39) =  0.000000000000000
-  M1066a_V%vs_1066a( 40) =  0.000000000000000
-  M1066a_V%vs_1066a( 41) =  0.000000000000000
-  M1066a_V%vs_1066a( 42) =  0.000000000000000
-  M1066a_V%vs_1066a( 43) =  0.000000000000000
-  M1066a_V%vs_1066a( 44) =  0.000000000000000
-  M1066a_V%vs_1066a( 45) =  0.000000000000000
-  M1066a_V%vs_1066a( 46) =  0.000000000000000
-  M1066a_V%vs_1066a( 47) =  0.000000000000000
-  M1066a_V%vs_1066a( 48) =  0.000000000000000
-  M1066a_V%vs_1066a( 49) =  0.000000000000000
-  M1066a_V%vs_1066a( 50) =  0.000000000000000
-  M1066a_V%vs_1066a( 51) =  0.000000000000000
-  M1066a_V%vs_1066a( 52) =  0.000000000000000
-  M1066a_V%vs_1066a( 53) =  0.000000000000000
-  M1066a_V%vs_1066a( 54) =  0.000000000000000
-  M1066a_V%vs_1066a( 55) =  0.000000000000000
-  M1066a_V%vs_1066a( 56) =  0.000000000000000
-  M1066a_V%vs_1066a( 57) =  0.000000000000000
-  M1066a_V%vs_1066a( 58) =  0.000000000000000
-  M1066a_V%vs_1066a( 59) =  0.000000000000000
-  M1066a_V%vs_1066a( 60) =  0.000000000000000
-  M1066a_V%vs_1066a( 61) =  0.000000000000000
-  M1066a_V%vs_1066a( 62) =  0.000000000000000
-  M1066a_V%vs_1066a( 63) =  0.000000000000000
-  M1066a_V%vs_1066a( 64) =  0.000000000000000
-  M1066a_V%vs_1066a( 65) =  0.000000000000000
-  M1066a_V%vs_1066a( 66) =  0.000000000000000
-  M1066a_V%vs_1066a( 67) =   7.24980000000000
-  M1066a_V%vs_1066a( 68) =   7.23760000000000
-  M1066a_V%vs_1066a( 69) =   7.22390000000000
-  M1066a_V%vs_1066a( 70) =   7.21000000000000
-  M1066a_V%vs_1066a( 71) =   7.19640000000000
-  M1066a_V%vs_1066a( 72) =   7.18300000000000
-  M1066a_V%vs_1066a( 73) =   7.16990000000000
-  M1066a_V%vs_1066a( 74) =   7.15710000000000
-  M1066a_V%vs_1066a( 75) =   7.14450000000000
-  M1066a_V%vs_1066a( 76) =   7.13200000000000
-  M1066a_V%vs_1066a( 77) =   7.11960000000000
-  M1066a_V%vs_1066a( 78) =   7.10740000000000
-  M1066a_V%vs_1066a( 79) =   7.09530000000000
-  M1066a_V%vs_1066a( 80) =   7.08320000000000
-  M1066a_V%vs_1066a( 81) =   7.07120000000000
-  M1066a_V%vs_1066a( 82) =   7.05920000000000
-  M1066a_V%vs_1066a( 83) =   7.04710000000000
-  M1066a_V%vs_1066a( 84) =   7.03470000000000
-  M1066a_V%vs_1066a( 85) =   7.02190000000000
-  M1066a_V%vs_1066a( 86) =   7.00860000000000
-  M1066a_V%vs_1066a( 87) =   6.99470000000000
-  M1066a_V%vs_1066a( 88) =   6.98030000000000
-  M1066a_V%vs_1066a( 89) =   6.96510000000000
-  M1066a_V%vs_1066a( 90) =   6.94930000000000
-  M1066a_V%vs_1066a( 91) =   6.93290000000000
-  M1066a_V%vs_1066a( 92) =   6.91620000000000
-  M1066a_V%vs_1066a( 93) =   6.89910000000000
-  M1066a_V%vs_1066a( 94) =   6.88200000000000
-  M1066a_V%vs_1066a( 95) =   6.86520000000000
-  M1066a_V%vs_1066a( 96) =   6.84900000000000
-  M1066a_V%vs_1066a( 97) =   6.83340000000000
-  M1066a_V%vs_1066a( 98) =   6.81820000000000
-  M1066a_V%vs_1066a( 99) =   6.80360000000000
-  M1066a_V%vs_1066a(100) =   6.78910000000000
-  M1066a_V%vs_1066a(101) =   6.77440000000000
-  M1066a_V%vs_1066a(102) =   6.75890000000000
-  M1066a_V%vs_1066a(103) =   6.74270000000000
-  M1066a_V%vs_1066a(104) =   6.72550000000000
-  M1066a_V%vs_1066a(105) =   6.70730000000000
-  M1066a_V%vs_1066a(106) =   6.68810000000000
-  M1066a_V%vs_1066a(107) =   6.66840000000000
-  M1066a_V%vs_1066a(108) =   6.64850000000000
-  M1066a_V%vs_1066a(109) =   6.62880000000000
-  M1066a_V%vs_1066a(110) =   6.60950000000000
-  M1066a_V%vs_1066a(111) =   6.59110000000000
-  M1066a_V%vs_1066a(112) =   6.57310000000000
-  M1066a_V%vs_1066a(113) =   6.55480000000000
-  M1066a_V%vs_1066a(114) =   6.53510000000000
-  M1066a_V%vs_1066a(115) =   6.51330000000000
-  M1066a_V%vs_1066a(116) =   6.48810000000000
-  M1066a_V%vs_1066a(117) =   6.45940000000000
-  M1066a_V%vs_1066a(118) =   6.42860000000000
-  M1066a_V%vs_1066a(119) =   6.39760000000000
-  M1066a_V%vs_1066a(120) =   6.36840000000000
-  M1066a_V%vs_1066a(121) =   6.34280000000000
-  M1066a_V%vs_1066a(122) =   6.32350000000000
-  M1066a_V%vs_1066a(123) =   6.31140000000000
-  M1066a_V%vs_1066a(124) =   6.30410000000000
-  M1066a_V%vs_1066a(125) =   6.30520000000000
-  M1066a_V%vs_1066a(126) =   6.30210000000000
-  M1066a_V%vs_1066a(127) =   6.26430000000000
-  M1066a_V%vs_1066a(128) =   6.19470000000000
-  M1066a_V%vs_1066a(129) =   6.09120000000000
-  M1066a_V%vs_1066a(130) =   5.95550000000000
-  M1066a_V%vs_1066a(131) =   5.77550000000000
-  M1066a_V%vs_1066a(132) =   5.77550000000000
-  M1066a_V%vs_1066a(133) =   5.60830000000000
-  M1066a_V%vs_1066a(134) =   5.47520000000000
-  M1066a_V%vs_1066a(135) =   5.36530000000000
-  M1066a_V%vs_1066a(136) =   5.26650000000000
-  M1066a_V%vs_1066a(137) =   5.17620000000000
-  M1066a_V%vs_1066a(138) =   5.09960000000000
-  M1066a_V%vs_1066a(139) =   5.03220000000000
-  M1066a_V%vs_1066a(140) =   4.94880000000000
-  M1066a_V%vs_1066a(141) =   4.94880000000000
-  M1066a_V%vs_1066a(142) =   4.86670000000000
-  M1066a_V%vs_1066a(143) =   4.78060000000000
-  M1066a_V%vs_1066a(144) =   4.69950000000000
-  M1066a_V%vs_1066a(145) =   4.62110000000000
-  M1066a_V%vs_1066a(146) =   4.54790000000000
-  M1066a_V%vs_1066a(147) =   4.48820000000000
-  M1066a_V%vs_1066a(148) =   4.44210000000000
-  M1066a_V%vs_1066a(149) =   4.40840000000000
-  M1066a_V%vs_1066a(150) =   4.38740000000000
-  M1066a_V%vs_1066a(151) =   4.37950000000000
-  M1066a_V%vs_1066a(152) =   4.39040000000000
-  M1066a_V%vs_1066a(153) =   4.43310000000000
-  M1066a_V%vs_1066a(154) =   4.48300000000000
-  M1066a_V%vs_1066a(155) =   4.53890000000000
-  M1066a_V%vs_1066a(156) =   4.60400000000000
-  M1066a_V%vs_1066a(157) =   4.64870000000000
-  M1066a_V%vs_1066a(158) =   2.58060000000000
-  M1066a_V%vs_1066a(159) =   2.58140000000000
-  M1066a_V%vs_1066a(160) =   2.58220000000000
+  M1066a_V%vs_1066a(  1) =   3.6298d0
+  M1066a_V%vs_1066a(  2) =   3.6297d0
+  M1066a_V%vs_1066a(  3) =   3.6294d0
+  M1066a_V%vs_1066a(  4) =   3.6288d0
+  M1066a_V%vs_1066a(  5) =   3.6281d0
+  M1066a_V%vs_1066a(  6) =   3.6271d0
+  M1066a_V%vs_1066a(  7) =   3.6259d0
+  M1066a_V%vs_1066a(  8) =   3.6244d0
+  M1066a_V%vs_1066a(  9) =   3.6228d0
+  M1066a_V%vs_1066a( 10) =   3.6209d0
+  M1066a_V%vs_1066a( 11) =   3.6187d0
+  M1066a_V%vs_1066a( 12) =   3.6163d0
+  M1066a_V%vs_1066a( 13) =   3.6137d0
+  M1066a_V%vs_1066a( 14) =   3.6108d0
+  M1066a_V%vs_1066a( 15) =   3.6076d0
+  M1066a_V%vs_1066a( 16) =   3.6042d0
+  M1066a_V%vs_1066a( 17) =   3.6004d0
+  M1066a_V%vs_1066a( 18) =   3.5965d0
+  M1066a_V%vs_1066a( 19) =   3.5922d0
+  M1066a_V%vs_1066a( 20) =   3.5876d0
+  M1066a_V%vs_1066a( 21) =   3.5828d0
+  M1066a_V%vs_1066a( 22) =   3.5777d0
+  M1066a_V%vs_1066a( 23) =   3.5724d0
+  M1066a_V%vs_1066a( 24) =   3.5668d0
+  M1066a_V%vs_1066a( 25) =   3.561d0
+  M1066a_V%vs_1066a( 26) =   3.5551d0
+  M1066a_V%vs_1066a( 27) =   3.549d0
+  M1066a_V%vs_1066a( 28) =   3.5428d0
+  M1066a_V%vs_1066a( 29) =   3.5365d0
+  M1066a_V%vs_1066a( 30) =   3.5301d0
+  M1066a_V%vs_1066a( 31) =   3.5238d0
+  M1066a_V%vs_1066a( 32) =   3.5172d0
+  M1066a_V%vs_1066a( 33) =   3.5118d0
+  M1066a_V%vs_1066a( 34) =  0.d0
+  M1066a_V%vs_1066a( 35) =  0.d0
+  M1066a_V%vs_1066a( 36) =  0.d0
+  M1066a_V%vs_1066a( 37) =  0.d0
+  M1066a_V%vs_1066a( 38) =  0.d0
+  M1066a_V%vs_1066a( 39) =  0.d0
+  M1066a_V%vs_1066a( 40) =  0.d0
+  M1066a_V%vs_1066a( 41) =  0.d0
+  M1066a_V%vs_1066a( 42) =  0.d0
+  M1066a_V%vs_1066a( 43) =  0.d0
+  M1066a_V%vs_1066a( 44) =  0.d0
+  M1066a_V%vs_1066a( 45) =  0.d0
+  M1066a_V%vs_1066a( 46) =  0.d0
+  M1066a_V%vs_1066a( 47) =  0.d0
+  M1066a_V%vs_1066a( 48) =  0.d0
+  M1066a_V%vs_1066a( 49) =  0.d0
+  M1066a_V%vs_1066a( 50) =  0.d0
+  M1066a_V%vs_1066a( 51) =  0.d0
+  M1066a_V%vs_1066a( 52) =  0.d0
+  M1066a_V%vs_1066a( 53) =  0.d0
+  M1066a_V%vs_1066a( 54) =  0.d0
+  M1066a_V%vs_1066a( 55) =  0.d0
+  M1066a_V%vs_1066a( 56) =  0.d0
+  M1066a_V%vs_1066a( 57) =  0.d0
+  M1066a_V%vs_1066a( 58) =  0.d0
+  M1066a_V%vs_1066a( 59) =  0.d0
+  M1066a_V%vs_1066a( 60) =  0.d0
+  M1066a_V%vs_1066a( 61) =  0.d0
+  M1066a_V%vs_1066a( 62) =  0.d0
+  M1066a_V%vs_1066a( 63) =  0.d0
+  M1066a_V%vs_1066a( 64) =  0.d0
+  M1066a_V%vs_1066a( 65) =  0.d0
+  M1066a_V%vs_1066a( 66) =  0.d0
+  M1066a_V%vs_1066a( 67) =   7.2498d0
+  M1066a_V%vs_1066a( 68) =   7.2376d0
+  M1066a_V%vs_1066a( 69) =   7.2239d0
+  M1066a_V%vs_1066a( 70) =   7.21d0
+  M1066a_V%vs_1066a( 71) =   7.1964d0
+  M1066a_V%vs_1066a( 72) =   7.183d0
+  M1066a_V%vs_1066a( 73) =   7.1699d0
+  M1066a_V%vs_1066a( 74) =   7.1571d0
+  M1066a_V%vs_1066a( 75) =   7.1445d0
+  M1066a_V%vs_1066a( 76) =   7.132d0
+  M1066a_V%vs_1066a( 77) =   7.1196d0
+  M1066a_V%vs_1066a( 78) =   7.1074d0
+  M1066a_V%vs_1066a( 79) =   7.0953d0
+  M1066a_V%vs_1066a( 80) =   7.0832d0
+  M1066a_V%vs_1066a( 81) =   7.0712d0
+  M1066a_V%vs_1066a( 82) =   7.0592d0
+  M1066a_V%vs_1066a( 83) =   7.0471d0
+  M1066a_V%vs_1066a( 84) =   7.0347d0
+  M1066a_V%vs_1066a( 85) =   7.0219d0
+  M1066a_V%vs_1066a( 86) =   7.0086d0
+  M1066a_V%vs_1066a( 87) =   6.9947d0
+  M1066a_V%vs_1066a( 88) =   6.9803d0
+  M1066a_V%vs_1066a( 89) =   6.9651d0
+  M1066a_V%vs_1066a( 90) =   6.9493d0
+  M1066a_V%vs_1066a( 91) =   6.9329d0
+  M1066a_V%vs_1066a( 92) =   6.9162d0
+  M1066a_V%vs_1066a( 93) =   6.8991d0
+  M1066a_V%vs_1066a( 94) =   6.882d0
+  M1066a_V%vs_1066a( 95) =   6.8652d0
+  M1066a_V%vs_1066a( 96) =   6.849d0
+  M1066a_V%vs_1066a( 97) =   6.8334d0
+  M1066a_V%vs_1066a( 98) =   6.8182d0
+  M1066a_V%vs_1066a( 99) =   6.8036d0
+  M1066a_V%vs_1066a(100) =   6.7891d0
+  M1066a_V%vs_1066a(101) =   6.7744d0
+  M1066a_V%vs_1066a(102) =   6.7589d0
+  M1066a_V%vs_1066a(103) =   6.7427d0
+  M1066a_V%vs_1066a(104) =   6.7255d0
+  M1066a_V%vs_1066a(105) =   6.7073d0
+  M1066a_V%vs_1066a(106) =   6.6881d0
+  M1066a_V%vs_1066a(107) =   6.6684d0
+  M1066a_V%vs_1066a(108) =   6.6485d0
+  M1066a_V%vs_1066a(109) =   6.6288d0
+  M1066a_V%vs_1066a(110) =   6.6095d0
+  M1066a_V%vs_1066a(111) =   6.5911d0
+  M1066a_V%vs_1066a(112) =   6.5731d0
+  M1066a_V%vs_1066a(113) =   6.5548d0
+  M1066a_V%vs_1066a(114) =   6.5351d0
+  M1066a_V%vs_1066a(115) =   6.5133d0
+  M1066a_V%vs_1066a(116) =   6.4881d0
+  M1066a_V%vs_1066a(117) =   6.4594d0
+  M1066a_V%vs_1066a(118) =   6.4286d0
+  M1066a_V%vs_1066a(119) =   6.3976d0
+  M1066a_V%vs_1066a(120) =   6.3684d0
+  M1066a_V%vs_1066a(121) =   6.3428d0
+  M1066a_V%vs_1066a(122) =   6.3235d0
+  M1066a_V%vs_1066a(123) =   6.3114d0
+  M1066a_V%vs_1066a(124) =   6.3041d0
+  M1066a_V%vs_1066a(125) =   6.3052d0
+  M1066a_V%vs_1066a(126) =   6.3021d0
+  M1066a_V%vs_1066a(127) =   6.2643d0
+  M1066a_V%vs_1066a(128) =   6.1947d0
+  M1066a_V%vs_1066a(129) =   6.0912d0
+  M1066a_V%vs_1066a(130) =   5.9555d0
+  M1066a_V%vs_1066a(131) =   5.7755d0
+  M1066a_V%vs_1066a(132) =   5.7755d0
+  M1066a_V%vs_1066a(133) =   5.6083d0
+  M1066a_V%vs_1066a(134) =   5.4752d0
+  M1066a_V%vs_1066a(135) =   5.3653d0
+  M1066a_V%vs_1066a(136) =   5.2665d0
+  M1066a_V%vs_1066a(137) =   5.1762d0
+  M1066a_V%vs_1066a(138) =   5.0996d0
+  M1066a_V%vs_1066a(139) =   5.0322d0
+  M1066a_V%vs_1066a(140) =   4.9488d0
+  M1066a_V%vs_1066a(141) =   4.9488d0
+  M1066a_V%vs_1066a(142) =   4.8667d0
+  M1066a_V%vs_1066a(143) =   4.7806d0
+  M1066a_V%vs_1066a(144) =   4.6995d0
+  M1066a_V%vs_1066a(145) =   4.6211d0
+  M1066a_V%vs_1066a(146) =   4.5479d0
+  M1066a_V%vs_1066a(147) =   4.4882d0
+  M1066a_V%vs_1066a(148) =   4.4421d0
+  M1066a_V%vs_1066a(149) =   4.4084d0
+  M1066a_V%vs_1066a(150) =   4.3874d0
+  M1066a_V%vs_1066a(151) =   4.3795d0
+  M1066a_V%vs_1066a(152) =   4.3904d0
+  M1066a_V%vs_1066a(153) =   4.4331d0
+  M1066a_V%vs_1066a(154) =   4.483d0
+  M1066a_V%vs_1066a(155) =   4.5389d0
+  M1066a_V%vs_1066a(156) =   4.604d0
+  M1066a_V%vs_1066a(157) =   4.6487d0
+  M1066a_V%vs_1066a(158) =   2.5806d0
+  M1066a_V%vs_1066a(159) =   2.5814d0
+  M1066a_V%vs_1066a(160) =   2.5822d0
 
   if (SUPPRESS_CRUSTAL_MESH) then
     M1066a_V%vp_1066a(158:160) = M1066a_V%vp_1066a(157)
@@ -835,327 +835,327 @@
     M1066a_V%density_1066a(158:160) = M1066a_V%density_1066a(157)
   endif
 
-  M1066a_V%Qkappa_1066a(  1) =   156900.000000000
-  M1066a_V%Qkappa_1066a(  2) =   156900.000000000
-  M1066a_V%Qkappa_1066a(  3) =   156900.000000000
-  M1066a_V%Qkappa_1066a(  4) =   156900.000000000
-  M1066a_V%Qkappa_1066a(  5) =   156900.000000000
-  M1066a_V%Qkappa_1066a(  6) =   156900.000000000
-  M1066a_V%Qkappa_1066a(  7) =   156900.000000000
-  M1066a_V%Qkappa_1066a(  8) =   156900.000000000
-  M1066a_V%Qkappa_1066a(  9) =   156900.000000000
-  M1066a_V%Qkappa_1066a( 10) =   156900.000000000
-  M1066a_V%Qkappa_1066a( 11) =   156900.000000000
-  M1066a_V%Qkappa_1066a( 12) =   156900.000000000
-  M1066a_V%Qkappa_1066a( 13) =   156900.000000000
-  M1066a_V%Qkappa_1066a( 14) =   156900.000000000
-  M1066a_V%Qkappa_1066a( 15) =   156900.000000000
-  M1066a_V%Qkappa_1066a( 16) =   156900.000000000
-  M1066a_V%Qkappa_1066a( 17) =   156900.000000000
-  M1066a_V%Qkappa_1066a( 18) =   156900.000000000
-  M1066a_V%Qkappa_1066a( 19) =   156900.000000000
-  M1066a_V%Qkappa_1066a( 20) =   156900.000000000
-  M1066a_V%Qkappa_1066a( 21) =   156900.000000000
-  M1066a_V%Qkappa_1066a( 22) =   156900.000000000
-  M1066a_V%Qkappa_1066a( 23) =   156900.000000000
-  M1066a_V%Qkappa_1066a( 24) =   156900.000000000
-  M1066a_V%Qkappa_1066a( 25) =   156900.000000000
-  M1066a_V%Qkappa_1066a( 26) =   156900.000000000
-  M1066a_V%Qkappa_1066a( 27) =   156900.000000000
-  M1066a_V%Qkappa_1066a( 28) =   156900.000000000
-  M1066a_V%Qkappa_1066a( 29) =   156900.000000000
-  M1066a_V%Qkappa_1066a( 30) =   156900.000000000
-  M1066a_V%Qkappa_1066a( 31) =   156900.000000000
-  M1066a_V%Qkappa_1066a( 32) =   156900.000000000
-  M1066a_V%Qkappa_1066a( 33) =   156900.000000000
-  M1066a_V%Qkappa_1066a( 34) =  0.000000000000000
-  M1066a_V%Qkappa_1066a( 35) =  0.000000000000000
-  M1066a_V%Qkappa_1066a( 36) =  0.000000000000000
-  M1066a_V%Qkappa_1066a( 37) =  0.000000000000000
-  M1066a_V%Qkappa_1066a( 38) =  0.000000000000000
-  M1066a_V%Qkappa_1066a( 39) =  0.000000000000000
-  M1066a_V%Qkappa_1066a( 40) =  0.000000000000000
-  M1066a_V%Qkappa_1066a( 41) =  0.000000000000000
-  M1066a_V%Qkappa_1066a( 42) =  0.000000000000000
-  M1066a_V%Qkappa_1066a( 43) =  0.000000000000000
-  M1066a_V%Qkappa_1066a( 44) =  0.000000000000000
-  M1066a_V%Qkappa_1066a( 45) =  0.000000000000000
-  M1066a_V%Qkappa_1066a( 46) =  0.000000000000000
-  M1066a_V%Qkappa_1066a( 47) =  0.000000000000000
-  M1066a_V%Qkappa_1066a( 48) =  0.000000000000000
-  M1066a_V%Qkappa_1066a( 49) =  0.000000000000000
-  M1066a_V%Qkappa_1066a( 50) =  0.000000000000000
-  M1066a_V%Qkappa_1066a( 51) =  0.000000000000000
-  M1066a_V%Qkappa_1066a( 52) =  0.000000000000000
-  M1066a_V%Qkappa_1066a( 53) =  0.000000000000000
-  M1066a_V%Qkappa_1066a( 54) =  0.000000000000000
-  M1066a_V%Qkappa_1066a( 55) =  0.000000000000000
-  M1066a_V%Qkappa_1066a( 56) =  0.000000000000000
-  M1066a_V%Qkappa_1066a( 57) =  0.000000000000000
-  M1066a_V%Qkappa_1066a( 58) =  0.000000000000000
-  M1066a_V%Qkappa_1066a( 59) =  0.000000000000000
-  M1066a_V%Qkappa_1066a( 60) =  0.000000000000000
-  M1066a_V%Qkappa_1066a( 61) =  0.000000000000000
-  M1066a_V%Qkappa_1066a( 62) =  0.000000000000000
-  M1066a_V%Qkappa_1066a( 63) =  0.000000000000000
-  M1066a_V%Qkappa_1066a( 64) =  0.000000000000000
-  M1066a_V%Qkappa_1066a( 65) =  0.000000000000000
-  M1066a_V%Qkappa_1066a( 66) =  0.000000000000000
-  M1066a_V%Qkappa_1066a( 67) =   16600.0000000000
-  M1066a_V%Qkappa_1066a( 68) =   16600.0000000000
-  M1066a_V%Qkappa_1066a( 69) =   16600.0000000000
-  M1066a_V%Qkappa_1066a( 70) =   16600.0000000000
-  M1066a_V%Qkappa_1066a( 71) =   16600.0000000000
-  M1066a_V%Qkappa_1066a( 72) =   16600.0000000000
-  M1066a_V%Qkappa_1066a( 73) =   16600.0000000000
-  M1066a_V%Qkappa_1066a( 74) =   16600.0000000000
-  M1066a_V%Qkappa_1066a( 75) =   16600.0000000000
-  M1066a_V%Qkappa_1066a( 76) =   16600.0000000000
-  M1066a_V%Qkappa_1066a( 77) =   16600.0000000000
-  M1066a_V%Qkappa_1066a( 78) =   16600.0000000000
-  M1066a_V%Qkappa_1066a( 79) =   16600.0000000000
-  M1066a_V%Qkappa_1066a( 80) =   16600.0000000000
-  M1066a_V%Qkappa_1066a( 81) =   16600.0000000000
-  M1066a_V%Qkappa_1066a( 82) =   16600.0000000000
-  M1066a_V%Qkappa_1066a( 83) =   16600.0000000000
-  M1066a_V%Qkappa_1066a( 84) =   16600.0000000000
-  M1066a_V%Qkappa_1066a( 85) =   16600.0000000000
-  M1066a_V%Qkappa_1066a( 86) =   16600.0000000000
-  M1066a_V%Qkappa_1066a( 87) =   16600.0000000000
-  M1066a_V%Qkappa_1066a( 88) =   16600.0000000000
-  M1066a_V%Qkappa_1066a( 89) =   16600.0000000000
-  M1066a_V%Qkappa_1066a( 90) =   16600.0000000000
-  M1066a_V%Qkappa_1066a( 91) =   16600.0000000000
-  M1066a_V%Qkappa_1066a( 92) =   16600.0000000000
-  M1066a_V%Qkappa_1066a( 93) =   16600.0000000000
-  M1066a_V%Qkappa_1066a( 94) =   16600.0000000000
-  M1066a_V%Qkappa_1066a( 95) =   16600.0000000000
-  M1066a_V%Qkappa_1066a( 96) =   16600.0000000000
-  M1066a_V%Qkappa_1066a( 97) =   16600.0000000000
-  M1066a_V%Qkappa_1066a( 98) =   16600.0000000000
-  M1066a_V%Qkappa_1066a( 99) =   16600.0000000000
-  M1066a_V%Qkappa_1066a(100) =   16600.0000000000
-  M1066a_V%Qkappa_1066a(101) =   16600.0000000000
-  M1066a_V%Qkappa_1066a(102) =   16600.0000000000
-  M1066a_V%Qkappa_1066a(103) =   16600.0000000000
-  M1066a_V%Qkappa_1066a(104) =   16600.0000000000
-  M1066a_V%Qkappa_1066a(105) =   16600.0000000000
-  M1066a_V%Qkappa_1066a(106) =   16600.0000000000
-  M1066a_V%Qkappa_1066a(107) =   16600.0000000000
-  M1066a_V%Qkappa_1066a(108) =   16600.0000000000
-  M1066a_V%Qkappa_1066a(109) =   16600.0000000000
-  M1066a_V%Qkappa_1066a(110) =   16600.0000000000
-  M1066a_V%Qkappa_1066a(111) =   16600.0000000000
-  M1066a_V%Qkappa_1066a(112) =   16600.0000000000
-  M1066a_V%Qkappa_1066a(113) =   16600.0000000000
-  M1066a_V%Qkappa_1066a(114) =   16600.0000000000
-  M1066a_V%Qkappa_1066a(115) =   16600.0000000000
-  M1066a_V%Qkappa_1066a(116) =   16600.0000000000
-  M1066a_V%Qkappa_1066a(117) =   16600.0000000000
-  M1066a_V%Qkappa_1066a(118) =   16600.0000000000
-  M1066a_V%Qkappa_1066a(119) =   16600.0000000000
-  M1066a_V%Qkappa_1066a(120) =   16600.0000000000
-  M1066a_V%Qkappa_1066a(121) =   16600.0000000000
-  M1066a_V%Qkappa_1066a(122) =   16600.0000000000
-  M1066a_V%Qkappa_1066a(123) =   16600.0000000000
-  M1066a_V%Qkappa_1066a(124) =   16600.0000000000
-  M1066a_V%Qkappa_1066a(125) =   16600.0000000000
-  M1066a_V%Qkappa_1066a(126) =   16600.0000000000
-  M1066a_V%Qkappa_1066a(127) =   16600.0000000000
-  M1066a_V%Qkappa_1066a(128) =   16600.0000000000
-  M1066a_V%Qkappa_1066a(129) =   16600.0000000000
-  M1066a_V%Qkappa_1066a(130) =   16600.0000000000
-  M1066a_V%Qkappa_1066a(131) =   16600.0000000000
-  M1066a_V%Qkappa_1066a(132) =   13840.0000000000
-  M1066a_V%Qkappa_1066a(133) =   13840.0000000000
-  M1066a_V%Qkappa_1066a(134) =   13840.0000000000
-  M1066a_V%Qkappa_1066a(135) =   13840.0000000000
-  M1066a_V%Qkappa_1066a(136) =   13840.0000000000
-  M1066a_V%Qkappa_1066a(137) =   13840.0000000000
-  M1066a_V%Qkappa_1066a(138) =   13840.0000000000
-  M1066a_V%Qkappa_1066a(139) =   13840.0000000000
-  M1066a_V%Qkappa_1066a(140) =   13840.0000000000
-  M1066a_V%Qkappa_1066a(141) =   5893.00000000000
-  M1066a_V%Qkappa_1066a(142) =   5893.00000000000
-  M1066a_V%Qkappa_1066a(143) =   5893.00000000000
-  M1066a_V%Qkappa_1066a(144) =   5893.00000000000
-  M1066a_V%Qkappa_1066a(145) =   5893.00000000000
-  M1066a_V%Qkappa_1066a(146) =   5893.00000000000
-  M1066a_V%Qkappa_1066a(147) =   5893.00000000000
-  M1066a_V%Qkappa_1066a(148) =   5893.00000000000
-  M1066a_V%Qkappa_1066a(149) =   5893.00000000000
-  M1066a_V%Qkappa_1066a(150) =   5893.00000000000
-  M1066a_V%Qkappa_1066a(151) =   5893.00000000000
-  M1066a_V%Qkappa_1066a(152) =   5893.00000000000
-  M1066a_V%Qkappa_1066a(153) =   5893.00000000000
-  M1066a_V%Qkappa_1066a(154) =   5893.00000000000
-  M1066a_V%Qkappa_1066a(155) =   5893.00000000000
-  M1066a_V%Qkappa_1066a(156) =   5893.00000000000
-  M1066a_V%Qkappa_1066a(157) =   5893.00000000000
-  M1066a_V%Qkappa_1066a(158) =   5893.00000000000
-  M1066a_V%Qkappa_1066a(159) =   5893.00000000000
-  M1066a_V%Qkappa_1066a(160) =   5893.00000000000
+  M1066a_V%Qkappa_1066a(  1) =   156900.d0
+  M1066a_V%Qkappa_1066a(  2) =   156900.d0
+  M1066a_V%Qkappa_1066a(  3) =   156900.d0
+  M1066a_V%Qkappa_1066a(  4) =   156900.d0
+  M1066a_V%Qkappa_1066a(  5) =   156900.d0
+  M1066a_V%Qkappa_1066a(  6) =   156900.d0
+  M1066a_V%Qkappa_1066a(  7) =   156900.d0
+  M1066a_V%Qkappa_1066a(  8) =   156900.d0
+  M1066a_V%Qkappa_1066a(  9) =   156900.d0
+  M1066a_V%Qkappa_1066a( 10) =   156900.d0
+  M1066a_V%Qkappa_1066a( 11) =   156900.d0
+  M1066a_V%Qkappa_1066a( 12) =   156900.d0
+  M1066a_V%Qkappa_1066a( 13) =   156900.d0
+  M1066a_V%Qkappa_1066a( 14) =   156900.d0
+  M1066a_V%Qkappa_1066a( 15) =   156900.d0
+  M1066a_V%Qkappa_1066a( 16) =   156900.d0
+  M1066a_V%Qkappa_1066a( 17) =   156900.d0
+  M1066a_V%Qkappa_1066a( 18) =   156900.d0
+  M1066a_V%Qkappa_1066a( 19) =   156900.d0
+  M1066a_V%Qkappa_1066a( 20) =   156900.d0
+  M1066a_V%Qkappa_1066a( 21) =   156900.d0
+  M1066a_V%Qkappa_1066a( 22) =   156900.d0
+  M1066a_V%Qkappa_1066a( 23) =   156900.d0
+  M1066a_V%Qkappa_1066a( 24) =   156900.d0
+  M1066a_V%Qkappa_1066a( 25) =   156900.d0
+  M1066a_V%Qkappa_1066a( 26) =   156900.d0
+  M1066a_V%Qkappa_1066a( 27) =   156900.d0
+  M1066a_V%Qkappa_1066a( 28) =   156900.d0
+  M1066a_V%Qkappa_1066a( 29) =   156900.d0
+  M1066a_V%Qkappa_1066a( 30) =   156900.d0
+  M1066a_V%Qkappa_1066a( 31) =   156900.d0
+  M1066a_V%Qkappa_1066a( 32) =   156900.d0
+  M1066a_V%Qkappa_1066a( 33) =   156900.d0
+  M1066a_V%Qkappa_1066a( 34) =  0.d0
+  M1066a_V%Qkappa_1066a( 35) =  0.d0
+  M1066a_V%Qkappa_1066a( 36) =  0.d0
+  M1066a_V%Qkappa_1066a( 37) =  0.d0
+  M1066a_V%Qkappa_1066a( 38) =  0.d0
+  M1066a_V%Qkappa_1066a( 39) =  0.d0
+  M1066a_V%Qkappa_1066a( 40) =  0.d0
+  M1066a_V%Qkappa_1066a( 41) =  0.d0
+  M1066a_V%Qkappa_1066a( 42) =  0.d0
+  M1066a_V%Qkappa_1066a( 43) =  0.d0
+  M1066a_V%Qkappa_1066a( 44) =  0.d0
+  M1066a_V%Qkappa_1066a( 45) =  0.d0
+  M1066a_V%Qkappa_1066a( 46) =  0.d0
+  M1066a_V%Qkappa_1066a( 47) =  0.d0
+  M1066a_V%Qkappa_1066a( 48) =  0.d0
+  M1066a_V%Qkappa_1066a( 49) =  0.d0
+  M1066a_V%Qkappa_1066a( 50) =  0.d0
+  M1066a_V%Qkappa_1066a( 51) =  0.d0
+  M1066a_V%Qkappa_1066a( 52) =  0.d0
+  M1066a_V%Qkappa_1066a( 53) =  0.d0
+  M1066a_V%Qkappa_1066a( 54) =  0.d0
+  M1066a_V%Qkappa_1066a( 55) =  0.d0
+  M1066a_V%Qkappa_1066a( 56) =  0.d0
+  M1066a_V%Qkappa_1066a( 57) =  0.d0
+  M1066a_V%Qkappa_1066a( 58) =  0.d0
+  M1066a_V%Qkappa_1066a( 59) =  0.d0
+  M1066a_V%Qkappa_1066a( 60) =  0.d0
+  M1066a_V%Qkappa_1066a( 61) =  0.d0
+  M1066a_V%Qkappa_1066a( 62) =  0.d0
+  M1066a_V%Qkappa_1066a( 63) =  0.d0
+  M1066a_V%Qkappa_1066a( 64) =  0.d0
+  M1066a_V%Qkappa_1066a( 65) =  0.d0
+  M1066a_V%Qkappa_1066a( 66) =  0.d0
+  M1066a_V%Qkappa_1066a( 67) =   16600.d0
+  M1066a_V%Qkappa_1066a( 68) =   16600.d0
+  M1066a_V%Qkappa_1066a( 69) =   16600.d0
+  M1066a_V%Qkappa_1066a( 70) =   16600.d0
+  M1066a_V%Qkappa_1066a( 71) =   16600.d0
+  M1066a_V%Qkappa_1066a( 72) =   16600.d0
+  M1066a_V%Qkappa_1066a( 73) =   16600.d0
+  M1066a_V%Qkappa_1066a( 74) =   16600.d0
+  M1066a_V%Qkappa_1066a( 75) =   16600.d0
+  M1066a_V%Qkappa_1066a( 76) =   16600.d0
+  M1066a_V%Qkappa_1066a( 77) =   16600.d0
+  M1066a_V%Qkappa_1066a( 78) =   16600.d0
+  M1066a_V%Qkappa_1066a( 79) =   16600.d0
+  M1066a_V%Qkappa_1066a( 80) =   16600.d0
+  M1066a_V%Qkappa_1066a( 81) =   16600.d0
+  M1066a_V%Qkappa_1066a( 82) =   16600.d0
+  M1066a_V%Qkappa_1066a( 83) =   16600.d0
+  M1066a_V%Qkappa_1066a( 84) =   16600.d0
+  M1066a_V%Qkappa_1066a( 85) =   16600.d0
+  M1066a_V%Qkappa_1066a( 86) =   16600.d0
+  M1066a_V%Qkappa_1066a( 87) =   16600.d0
+  M1066a_V%Qkappa_1066a( 88) =   16600.d0
+  M1066a_V%Qkappa_1066a( 89) =   16600.d0
+  M1066a_V%Qkappa_1066a( 90) =   16600.d0
+  M1066a_V%Qkappa_1066a( 91) =   16600.d0
+  M1066a_V%Qkappa_1066a( 92) =   16600.d0
+  M1066a_V%Qkappa_1066a( 93) =   16600.d0
+  M1066a_V%Qkappa_1066a( 94) =   16600.d0
+  M1066a_V%Qkappa_1066a( 95) =   16600.d0
+  M1066a_V%Qkappa_1066a( 96) =   16600.d0
+  M1066a_V%Qkappa_1066a( 97) =   16600.d0
+  M1066a_V%Qkappa_1066a( 98) =   16600.d0
+  M1066a_V%Qkappa_1066a( 99) =   16600.d0
+  M1066a_V%Qkappa_1066a(100) =   16600.d0
+  M1066a_V%Qkappa_1066a(101) =   16600.d0
+  M1066a_V%Qkappa_1066a(102) =   16600.d0
+  M1066a_V%Qkappa_1066a(103) =   16600.d0
+  M1066a_V%Qkappa_1066a(104) =   16600.d0
+  M1066a_V%Qkappa_1066a(105) =   16600.d0
+  M1066a_V%Qkappa_1066a(106) =   16600.d0
+  M1066a_V%Qkappa_1066a(107) =   16600.d0
+  M1066a_V%Qkappa_1066a(108) =   16600.d0
+  M1066a_V%Qkappa_1066a(109) =   16600.d0
+  M1066a_V%Qkappa_1066a(110) =   16600.d0
+  M1066a_V%Qkappa_1066a(111) =   16600.d0
+  M1066a_V%Qkappa_1066a(112) =   16600.d0
+  M1066a_V%Qkappa_1066a(113) =   16600.d0
+  M1066a_V%Qkappa_1066a(114) =   16600.d0
+  M1066a_V%Qkappa_1066a(115) =   16600.d0
+  M1066a_V%Qkappa_1066a(116) =   16600.d0
+  M1066a_V%Qkappa_1066a(117) =   16600.d0
+  M1066a_V%Qkappa_1066a(118) =   16600.d0
+  M1066a_V%Qkappa_1066a(119) =   16600.d0
+  M1066a_V%Qkappa_1066a(120) =   16600.d0
+  M1066a_V%Qkappa_1066a(121) =   16600.d0
+  M1066a_V%Qkappa_1066a(122) =   16600.d0
+  M1066a_V%Qkappa_1066a(123) =   16600.d0
+  M1066a_V%Qkappa_1066a(124) =   16600.d0
+  M1066a_V%Qkappa_1066a(125) =   16600.d0
+  M1066a_V%Qkappa_1066a(126) =   16600.d0
+  M1066a_V%Qkappa_1066a(127) =   16600.d0
+  M1066a_V%Qkappa_1066a(128) =   16600.d0
+  M1066a_V%Qkappa_1066a(129) =   16600.d0
+  M1066a_V%Qkappa_1066a(130) =   16600.d0
+  M1066a_V%Qkappa_1066a(131) =   16600.d0
+  M1066a_V%Qkappa_1066a(132) =   13840.d0
+  M1066a_V%Qkappa_1066a(133) =   13840.d0
+  M1066a_V%Qkappa_1066a(134) =   13840.d0
+  M1066a_V%Qkappa_1066a(135) =   13840.d0
+  M1066a_V%Qkappa_1066a(136) =   13840.d0
+  M1066a_V%Qkappa_1066a(137) =   13840.d0
+  M1066a_V%Qkappa_1066a(138) =   13840.d0
+  M1066a_V%Qkappa_1066a(139) =   13840.d0
+  M1066a_V%Qkappa_1066a(140) =   13840.d0
+  M1066a_V%Qkappa_1066a(141) =   5893.d0
+  M1066a_V%Qkappa_1066a(142) =   5893.d0
+  M1066a_V%Qkappa_1066a(143) =   5893.d0
+  M1066a_V%Qkappa_1066a(144) =   5893.d0
+  M1066a_V%Qkappa_1066a(145) =   5893.d0
+  M1066a_V%Qkappa_1066a(146) =   5893.d0
+  M1066a_V%Qkappa_1066a(147) =   5893.d0
+  M1066a_V%Qkappa_1066a(148) =   5893.d0
+  M1066a_V%Qkappa_1066a(149) =   5893.d0
+  M1066a_V%Qkappa_1066a(150) =   5893.d0
+  M1066a_V%Qkappa_1066a(151) =   5893.d0
+  M1066a_V%Qkappa_1066a(152) =   5893.d0
+  M1066a_V%Qkappa_1066a(153) =   5893.d0
+  M1066a_V%Qkappa_1066a(154) =   5893.d0
+  M1066a_V%Qkappa_1066a(155) =   5893.d0
+  M1066a_V%Qkappa_1066a(156) =   5893.d0
+  M1066a_V%Qkappa_1066a(157) =   5893.d0
+  M1066a_V%Qkappa_1066a(158) =   5893.d0
+  M1066a_V%Qkappa_1066a(159) =   5893.d0
+  M1066a_V%Qkappa_1066a(160) =   5893.d0
 
-  M1066a_V%Qmu_1066a(  1) =   3138.00000000000
-  M1066a_V%Qmu_1066a(  2) =   3138.00000000000
-  M1066a_V%Qmu_1066a(  3) =   3138.00000000000
-  M1066a_V%Qmu_1066a(  4) =   3138.00000000000
-  M1066a_V%Qmu_1066a(  5) =   3138.00000000000
-  M1066a_V%Qmu_1066a(  6) =   3138.00000000000
-  M1066a_V%Qmu_1066a(  7) =   3138.00000000000
-  M1066a_V%Qmu_1066a(  8) =   3138.00000000000
-  M1066a_V%Qmu_1066a(  9) =   3138.00000000000
-  M1066a_V%Qmu_1066a( 10) =   3138.00000000000
-  M1066a_V%Qmu_1066a( 11) =   3138.00000000000
-  M1066a_V%Qmu_1066a( 12) =   3138.00000000000
-  M1066a_V%Qmu_1066a( 13) =   3138.00000000000
-  M1066a_V%Qmu_1066a( 14) =   3138.00000000000
-  M1066a_V%Qmu_1066a( 15) =   3138.00000000000
-  M1066a_V%Qmu_1066a( 16) =   3138.00000000000
-  M1066a_V%Qmu_1066a( 17) =   3138.00000000000
-  M1066a_V%Qmu_1066a( 18) =   3138.00000000000
-  M1066a_V%Qmu_1066a( 19) =   3138.00000000000
-  M1066a_V%Qmu_1066a( 20) =   3138.00000000000
-  M1066a_V%Qmu_1066a( 21) =   3138.00000000000
-  M1066a_V%Qmu_1066a( 22) =   3138.00000000000
-  M1066a_V%Qmu_1066a( 23) =   3138.00000000000
-  M1066a_V%Qmu_1066a( 24) =   3138.00000000000
-  M1066a_V%Qmu_1066a( 25) =   3138.00000000000
-  M1066a_V%Qmu_1066a( 26) =   3138.00000000000
-  M1066a_V%Qmu_1066a( 27) =   3138.00000000000
-  M1066a_V%Qmu_1066a( 28) =   3138.00000000000
-  M1066a_V%Qmu_1066a( 29) =   3138.00000000000
-  M1066a_V%Qmu_1066a( 30) =   3138.00000000000
-  M1066a_V%Qmu_1066a( 31) =   3138.00000000000
-  M1066a_V%Qmu_1066a( 32) =   3138.00000000000
-  M1066a_V%Qmu_1066a( 33) =   3138.00000000000
-  M1066a_V%Qmu_1066a( 34) =  0.000000000000000
-  M1066a_V%Qmu_1066a( 35) =  0.000000000000000
-  M1066a_V%Qmu_1066a( 36) =  0.000000000000000
-  M1066a_V%Qmu_1066a( 37) =  0.000000000000000
-  M1066a_V%Qmu_1066a( 38) =  0.000000000000000
-  M1066a_V%Qmu_1066a( 39) =  0.000000000000000
-  M1066a_V%Qmu_1066a( 40) =  0.000000000000000
-  M1066a_V%Qmu_1066a( 41) =  0.000000000000000
-  M1066a_V%Qmu_1066a( 42) =  0.000000000000000
-  M1066a_V%Qmu_1066a( 43) =  0.000000000000000
-  M1066a_V%Qmu_1066a( 44) =  0.000000000000000
-  M1066a_V%Qmu_1066a( 45) =  0.000000000000000
-  M1066a_V%Qmu_1066a( 46) =  0.000000000000000
-  M1066a_V%Qmu_1066a( 47) =  0.000000000000000
-  M1066a_V%Qmu_1066a( 48) =  0.000000000000000
-  M1066a_V%Qmu_1066a( 49) =  0.000000000000000
-  M1066a_V%Qmu_1066a( 50) =  0.000000000000000
-  M1066a_V%Qmu_1066a( 51) =  0.000000000000000
-  M1066a_V%Qmu_1066a( 52) =  0.000000000000000
-  M1066a_V%Qmu_1066a( 53) =  0.000000000000000
-  M1066a_V%Qmu_1066a( 54) =  0.000000000000000
-  M1066a_V%Qmu_1066a( 55) =  0.000000000000000
-  M1066a_V%Qmu_1066a( 56) =  0.000000000000000
-  M1066a_V%Qmu_1066a( 57) =  0.000000000000000
-  M1066a_V%Qmu_1066a( 58) =  0.000000000000000
-  M1066a_V%Qmu_1066a( 59) =  0.000000000000000
-  M1066a_V%Qmu_1066a( 60) =  0.000000000000000
-  M1066a_V%Qmu_1066a( 61) =  0.000000000000000
-  M1066a_V%Qmu_1066a( 62) =  0.000000000000000
-  M1066a_V%Qmu_1066a( 63) =  0.000000000000000
-  M1066a_V%Qmu_1066a( 64) =  0.000000000000000
-  M1066a_V%Qmu_1066a( 65) =  0.000000000000000
-  M1066a_V%Qmu_1066a( 66) =  0.000000000000000
-  M1066a_V%Qmu_1066a( 67) =   332.000000000000
-  M1066a_V%Qmu_1066a( 68) =   332.000000000000
-  M1066a_V%Qmu_1066a( 69) =   332.000000000000
-  M1066a_V%Qmu_1066a( 70) =   332.000000000000
-  M1066a_V%Qmu_1066a( 71) =   332.000000000000
-  M1066a_V%Qmu_1066a( 72) =   332.000000000000
-  M1066a_V%Qmu_1066a( 73) =   332.000000000000
-  M1066a_V%Qmu_1066a( 74) =   332.000000000000
-  M1066a_V%Qmu_1066a( 75) =   332.000000000000
-  M1066a_V%Qmu_1066a( 76) =   332.000000000000
-  M1066a_V%Qmu_1066a( 77) =   332.000000000000
-  M1066a_V%Qmu_1066a( 78) =   332.000000000000
-  M1066a_V%Qmu_1066a( 79) =   332.000000000000
-  M1066a_V%Qmu_1066a( 80) =   332.000000000000
-  M1066a_V%Qmu_1066a( 81) =   332.000000000000
-  M1066a_V%Qmu_1066a( 82) =   332.000000000000
-  M1066a_V%Qmu_1066a( 83) =   332.000000000000
-  M1066a_V%Qmu_1066a( 84) =   332.000000000000
-  M1066a_V%Qmu_1066a( 85) =   332.000000000000
-  M1066a_V%Qmu_1066a( 86) =   332.000000000000
-  M1066a_V%Qmu_1066a( 87) =   332.000000000000
-  M1066a_V%Qmu_1066a( 88) =   332.000000000000
-  M1066a_V%Qmu_1066a( 89) =   332.000000000000
-  M1066a_V%Qmu_1066a( 90) =   332.000000000000
-  M1066a_V%Qmu_1066a( 91) =   332.000000000000
-  M1066a_V%Qmu_1066a( 92) =   332.000000000000
-  M1066a_V%Qmu_1066a( 93) =   332.000000000000
-  M1066a_V%Qmu_1066a( 94) =   332.000000000000
-  M1066a_V%Qmu_1066a( 95) =   332.000000000000
-  M1066a_V%Qmu_1066a( 96) =   332.000000000000
-  M1066a_V%Qmu_1066a( 97) =   332.000000000000
-  M1066a_V%Qmu_1066a( 98) =   332.000000000000
-  M1066a_V%Qmu_1066a( 99) =   332.000000000000
-  M1066a_V%Qmu_1066a(100) =   332.000000000000
-  M1066a_V%Qmu_1066a(101) =   332.000000000000
-  M1066a_V%Qmu_1066a(102) =   332.000000000000
-  M1066a_V%Qmu_1066a(103) =   332.000000000000
-  M1066a_V%Qmu_1066a(104) =   332.000000000000
-  M1066a_V%Qmu_1066a(105) =   332.000000000000
-  M1066a_V%Qmu_1066a(106) =   332.000000000000
-  M1066a_V%Qmu_1066a(107) =   332.000000000000
-  M1066a_V%Qmu_1066a(108) =   332.000000000000
-  M1066a_V%Qmu_1066a(109) =   332.000000000000
-  M1066a_V%Qmu_1066a(110) =   332.000000000000
-  M1066a_V%Qmu_1066a(111) =   332.000000000000
-  M1066a_V%Qmu_1066a(112) =   332.000000000000
-  M1066a_V%Qmu_1066a(113) =   332.000000000000
-  M1066a_V%Qmu_1066a(114) =   332.000000000000
-  M1066a_V%Qmu_1066a(115) =   332.000000000000
-  M1066a_V%Qmu_1066a(116) =   332.000000000000
-  M1066a_V%Qmu_1066a(117) =   332.000000000000
-  M1066a_V%Qmu_1066a(118) =   332.000000000000
-  M1066a_V%Qmu_1066a(119) =   332.000000000000
-  M1066a_V%Qmu_1066a(120) =   332.000000000000
-  M1066a_V%Qmu_1066a(121) =   332.000000000000
-  M1066a_V%Qmu_1066a(122) =   332.000000000000
-  M1066a_V%Qmu_1066a(123) =   332.000000000000
-  M1066a_V%Qmu_1066a(124) =   332.000000000000
-  M1066a_V%Qmu_1066a(125) =   332.000000000000
-  M1066a_V%Qmu_1066a(126) =   332.000000000000
-  M1066a_V%Qmu_1066a(127) =   332.000000000000
-  M1066a_V%Qmu_1066a(128) =   332.000000000000
-  M1066a_V%Qmu_1066a(129) =   332.000000000000
-  M1066a_V%Qmu_1066a(130) =   332.000000000000
-  M1066a_V%Qmu_1066a(131) =   332.000000000000
-  M1066a_V%Qmu_1066a(132) =   276.800000000000
-  M1066a_V%Qmu_1066a(133) =   276.800000000000
-  M1066a_V%Qmu_1066a(134) =   276.800000000000
-  M1066a_V%Qmu_1066a(135) =   276.800000000000
-  M1066a_V%Qmu_1066a(136) =   276.800000000000
-  M1066a_V%Qmu_1066a(137) =   276.800000000000
-  M1066a_V%Qmu_1066a(138) =   276.800000000000
-  M1066a_V%Qmu_1066a(139) =   276.800000000000
-  M1066a_V%Qmu_1066a(140) =   276.800000000000
-  M1066a_V%Qmu_1066a(141) =   117.900000000000
-  M1066a_V%Qmu_1066a(142) =   117.900000000000
-  M1066a_V%Qmu_1066a(143) =   117.900000000000
-  M1066a_V%Qmu_1066a(144) =   117.900000000000
-  M1066a_V%Qmu_1066a(145) =   117.900000000000
-  M1066a_V%Qmu_1066a(146) =   117.900000000000
-  M1066a_V%Qmu_1066a(147) =   117.900000000000
-  M1066a_V%Qmu_1066a(148) =   117.900000000000
-  M1066a_V%Qmu_1066a(149) =   117.900000000000
-  M1066a_V%Qmu_1066a(150) =   117.900000000000
-  M1066a_V%Qmu_1066a(151) =   117.900000000000
-  M1066a_V%Qmu_1066a(152) =   117.900000000000
-  M1066a_V%Qmu_1066a(153) =   117.900000000000
-  M1066a_V%Qmu_1066a(154) =   117.900000000000
-  M1066a_V%Qmu_1066a(155) =   117.900000000000
-  M1066a_V%Qmu_1066a(156) =   117.900000000000
-  M1066a_V%Qmu_1066a(157) =   117.900000000000
-  M1066a_V%Qmu_1066a(158) =   117.900000000000
-  M1066a_V%Qmu_1066a(159) =   117.900000000000
-  M1066a_V%Qmu_1066a(160) =   117.900000000000
+  M1066a_V%Qmu_1066a(  1) =   3138.d0
+  M1066a_V%Qmu_1066a(  2) =   3138.d0
+  M1066a_V%Qmu_1066a(  3) =   3138.d0
+  M1066a_V%Qmu_1066a(  4) =   3138.d0
+  M1066a_V%Qmu_1066a(  5) =   3138.d0
+  M1066a_V%Qmu_1066a(  6) =   3138.d0
+  M1066a_V%Qmu_1066a(  7) =   3138.d0
+  M1066a_V%Qmu_1066a(  8) =   3138.d0
+  M1066a_V%Qmu_1066a(  9) =   3138.d0
+  M1066a_V%Qmu_1066a( 10) =   3138.d0
+  M1066a_V%Qmu_1066a( 11) =   3138.d0
+  M1066a_V%Qmu_1066a( 12) =   3138.d0
+  M1066a_V%Qmu_1066a( 13) =   3138.d0
+  M1066a_V%Qmu_1066a( 14) =   3138.d0
+  M1066a_V%Qmu_1066a( 15) =   3138.d0
+  M1066a_V%Qmu_1066a( 16) =   3138.d0
+  M1066a_V%Qmu_1066a( 17) =   3138.d0
+  M1066a_V%Qmu_1066a( 18) =   3138.d0
+  M1066a_V%Qmu_1066a( 19) =   3138.d0
+  M1066a_V%Qmu_1066a( 20) =   3138.d0
+  M1066a_V%Qmu_1066a( 21) =   3138.d0
+  M1066a_V%Qmu_1066a( 22) =   3138.d0
+  M1066a_V%Qmu_1066a( 23) =   3138.d0
+  M1066a_V%Qmu_1066a( 24) =   3138.d0
+  M1066a_V%Qmu_1066a( 25) =   3138.d0
+  M1066a_V%Qmu_1066a( 26) =   3138.d0
+  M1066a_V%Qmu_1066a( 27) =   3138.d0
+  M1066a_V%Qmu_1066a( 28) =   3138.d0
+  M1066a_V%Qmu_1066a( 29) =   3138.d0
+  M1066a_V%Qmu_1066a( 30) =   3138.d0
+  M1066a_V%Qmu_1066a( 31) =   3138.d0
+  M1066a_V%Qmu_1066a( 32) =   3138.d0
+  M1066a_V%Qmu_1066a( 33) =   3138.d0
+  M1066a_V%Qmu_1066a( 34) =  0.d0
+  M1066a_V%Qmu_1066a( 35) =  0.d0
+  M1066a_V%Qmu_1066a( 36) =  0.d0
+  M1066a_V%Qmu_1066a( 37) =  0.d0
+  M1066a_V%Qmu_1066a( 38) =  0.d0
+  M1066a_V%Qmu_1066a( 39) =  0.d0
+  M1066a_V%Qmu_1066a( 40) =  0.d0
+  M1066a_V%Qmu_1066a( 41) =  0.d0
+  M1066a_V%Qmu_1066a( 42) =  0.d0
+  M1066a_V%Qmu_1066a( 43) =  0.d0
+  M1066a_V%Qmu_1066a( 44) =  0.d0
+  M1066a_V%Qmu_1066a( 45) =  0.d0
+  M1066a_V%Qmu_1066a( 46) =  0.d0
+  M1066a_V%Qmu_1066a( 47) =  0.d0
+  M1066a_V%Qmu_1066a( 48) =  0.d0
+  M1066a_V%Qmu_1066a( 49) =  0.d0
+  M1066a_V%Qmu_1066a( 50) =  0.d0
+  M1066a_V%Qmu_1066a( 51) =  0.d0
+  M1066a_V%Qmu_1066a( 52) =  0.d0
+  M1066a_V%Qmu_1066a( 53) =  0.d0
+  M1066a_V%Qmu_1066a( 54) =  0.d0
+  M1066a_V%Qmu_1066a( 55) =  0.d0
+  M1066a_V%Qmu_1066a( 56) =  0.d0
+  M1066a_V%Qmu_1066a( 57) =  0.d0
+  M1066a_V%Qmu_1066a( 58) =  0.d0
+  M1066a_V%Qmu_1066a( 59) =  0.d0
+  M1066a_V%Qmu_1066a( 60) =  0.d0
+  M1066a_V%Qmu_1066a( 61) =  0.d0
+  M1066a_V%Qmu_1066a( 62) =  0.d0
+  M1066a_V%Qmu_1066a( 63) =  0.d0
+  M1066a_V%Qmu_1066a( 64) =  0.d0
+  M1066a_V%Qmu_1066a( 65) =  0.d0
+  M1066a_V%Qmu_1066a( 66) =  0.d0
+  M1066a_V%Qmu_1066a( 67) =   332.d0
+  M1066a_V%Qmu_1066a( 68) =   332.d0
+  M1066a_V%Qmu_1066a( 69) =   332.d0
+  M1066a_V%Qmu_1066a( 70) =   332.d0
+  M1066a_V%Qmu_1066a( 71) =   332.d0
+  M1066a_V%Qmu_1066a( 72) =   332.d0
+  M1066a_V%Qmu_1066a( 73) =   332.d0
+  M1066a_V%Qmu_1066a( 74) =   332.d0
+  M1066a_V%Qmu_1066a( 75) =   332.d0
+  M1066a_V%Qmu_1066a( 76) =   332.d0
+  M1066a_V%Qmu_1066a( 77) =   332.d0
+  M1066a_V%Qmu_1066a( 78) =   332.d0
+  M1066a_V%Qmu_1066a( 79) =   332.d0
+  M1066a_V%Qmu_1066a( 80) =   332.d0
+  M1066a_V%Qmu_1066a( 81) =   332.d0
+  M1066a_V%Qmu_1066a( 82) =   332.d0
+  M1066a_V%Qmu_1066a( 83) =   332.d0
+  M1066a_V%Qmu_1066a( 84) =   332.d0
+  M1066a_V%Qmu_1066a( 85) =   332.d0
+  M1066a_V%Qmu_1066a( 86) =   332.d0
+  M1066a_V%Qmu_1066a( 87) =   332.d0
+  M1066a_V%Qmu_1066a( 88) =   332.d0
+  M1066a_V%Qmu_1066a( 89) =   332.d0
+  M1066a_V%Qmu_1066a( 90) =   332.d0
+  M1066a_V%Qmu_1066a( 91) =   332.d0
+  M1066a_V%Qmu_1066a( 92) =   332.d0
+  M1066a_V%Qmu_1066a( 93) =   332.d0
+  M1066a_V%Qmu_1066a( 94) =   332.d0
+  M1066a_V%Qmu_1066a( 95) =   332.d0
+  M1066a_V%Qmu_1066a( 96) =   332.d0
+  M1066a_V%Qmu_1066a( 97) =   332.d0
+  M1066a_V%Qmu_1066a( 98) =   332.d0
+  M1066a_V%Qmu_1066a( 99) =   332.d0
+  M1066a_V%Qmu_1066a(100) =   332.d0
+  M1066a_V%Qmu_1066a(101) =   332.d0
+  M1066a_V%Qmu_1066a(102) =   332.d0
+  M1066a_V%Qmu_1066a(103) =   332.d0
+  M1066a_V%Qmu_1066a(104) =   332.d0
+  M1066a_V%Qmu_1066a(105) =   332.d0
+  M1066a_V%Qmu_1066a(106) =   332.d0
+  M1066a_V%Qmu_1066a(107) =   332.d0
+  M1066a_V%Qmu_1066a(108) =   332.d0
+  M1066a_V%Qmu_1066a(109) =   332.d0
+  M1066a_V%Qmu_1066a(110) =   332.d0
+  M1066a_V%Qmu_1066a(111) =   332.d0
+  M1066a_V%Qmu_1066a(112) =   332.d0
+  M1066a_V%Qmu_1066a(113) =   332.d0
+  M1066a_V%Qmu_1066a(114) =   332.d0
+  M1066a_V%Qmu_1066a(115) =   332.d0
+  M1066a_V%Qmu_1066a(116) =   332.d0
+  M1066a_V%Qmu_1066a(117) =   332.d0
+  M1066a_V%Qmu_1066a(118) =   332.d0
+  M1066a_V%Qmu_1066a(119) =   332.d0
+  M1066a_V%Qmu_1066a(120) =   332.d0
+  M1066a_V%Qmu_1066a(121) =   332.d0
+  M1066a_V%Qmu_1066a(122) =   332.d0
+  M1066a_V%Qmu_1066a(123) =   332.d0
+  M1066a_V%Qmu_1066a(124) =   332.d0
+  M1066a_V%Qmu_1066a(125) =   332.d0
+  M1066a_V%Qmu_1066a(126) =   332.d0
+  M1066a_V%Qmu_1066a(127) =   332.d0
+  M1066a_V%Qmu_1066a(128) =   332.d0
+  M1066a_V%Qmu_1066a(129) =   332.d0
+  M1066a_V%Qmu_1066a(130) =   332.d0
+  M1066a_V%Qmu_1066a(131) =   332.d0
+  M1066a_V%Qmu_1066a(132) =   276.8d0
+  M1066a_V%Qmu_1066a(133) =   276.8d0
+  M1066a_V%Qmu_1066a(134) =   276.8d0
+  M1066a_V%Qmu_1066a(135) =   276.8d0
+  M1066a_V%Qmu_1066a(136) =   276.8d0
+  M1066a_V%Qmu_1066a(137) =   276.8d0
+  M1066a_V%Qmu_1066a(138) =   276.8d0
+  M1066a_V%Qmu_1066a(139) =   276.8d0
+  M1066a_V%Qmu_1066a(140) =   276.8d0
+  M1066a_V%Qmu_1066a(141) =   117.9d0
+  M1066a_V%Qmu_1066a(142) =   117.9d0
+  M1066a_V%Qmu_1066a(143) =   117.9d0
+  M1066a_V%Qmu_1066a(144) =   117.9d0
+  M1066a_V%Qmu_1066a(145) =   117.9d0
+  M1066a_V%Qmu_1066a(146) =   117.9d0
+  M1066a_V%Qmu_1066a(147) =   117.9d0
+  M1066a_V%Qmu_1066a(148) =   117.9d0
+  M1066a_V%Qmu_1066a(149) =   117.9d0
+  M1066a_V%Qmu_1066a(150) =   117.9d0
+  M1066a_V%Qmu_1066a(151) =   117.9d0
+  M1066a_V%Qmu_1066a(152) =   117.9d0
+  M1066a_V%Qmu_1066a(153) =   117.9d0
+  M1066a_V%Qmu_1066a(154) =   117.9d0
+  M1066a_V%Qmu_1066a(155) =   117.9d0
+  M1066a_V%Qmu_1066a(156) =   117.9d0
+  M1066a_V%Qmu_1066a(157) =   117.9d0
+  M1066a_V%Qmu_1066a(158) =   117.9d0
+  M1066a_V%Qmu_1066a(159) =   117.9d0
+  M1066a_V%Qmu_1066a(160) =   117.9d0
 
 ! strip the crust and replace it by mantle if we use an external crustal model
   if (SUPPRESS_CRUSTAL_MESH .or. USE_EXTERNAL_CRUSTAL_MODEL) then

Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/model_1dref.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/model_1dref.f90	2013-07-24 14:30:04 UTC (rev 22661)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/model_1dref.f90	2013-07-24 15:33:12 UTC (rev 22662)
@@ -223,7 +223,7 @@
   logical USE_EXTERNAL_CRUSTAL_MODEL
 
 
-! define the 1D REF model of Kustowski et al. (2007)
+  ! define the 1D REF model of Kustowski et al. (2007)
 
  Mref_V%radius_ref( 1 : 30 ) = (/ &
  0.000000000000000E+000 , &

Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/model_ak135.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/model_ak135.f90	2013-07-24 14:30:04 UTC (rev 22661)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/model_ak135.f90	2013-07-24 15:33:12 UTC (rev 22662)
@@ -195,7 +195,7 @@
  type (model_ak135_variables) Mak135_V
 ! model_ak135_variables
 
-  logical USE_EXTERNAL_CRUSTAL_MODEL
+  logical :: USE_EXTERNAL_CRUSTAL_MODEL
 
 ! define all the values in the model
 

Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/model_aniso_mantle.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/model_aniso_mantle.f90	2013-07-24 14:30:04 UTC (rev 22661)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/model_aniso_mantle.f90	2013-07-24 15:33:12 UTC (rev 22662)
@@ -64,6 +64,8 @@
   ! model_aniso_mantle_variables
 
   integer :: myrank
+
+  ! local parameters
   integer :: ier
 
   ! the variables read are declared and stored in structure AMM_V

Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/model_atten3D_QRFSI12.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/model_atten3D_QRFSI12.f90	2013-07-24 14:30:04 UTC (rev 22661)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/model_atten3D_QRFSI12.f90	2013-07-24 15:33:12 UTC (rev 22662)
@@ -149,14 +149,14 @@
   close(10)
 
 ! get the depths (km) of the spline knots
-  QRFSI12_Q%spknt(1) = 24.4
-  QRFSI12_Q%spknt(2) = 75.0
-  QRFSI12_Q%spknt(3) = 150.0
-  QRFSI12_Q%spknt(4) = 225.0
-  QRFSI12_Q%spknt(5) = 300.0
-  QRFSI12_Q%spknt(6) = 410.0
-  QRFSI12_Q%spknt(7) = 530.0
-  QRFSI12_Q%spknt(8) = 650.0
+  QRFSI12_Q%spknt(1) = 24.4d0
+  QRFSI12_Q%spknt(2) = 75.d0
+  QRFSI12_Q%spknt(3) = 150.d0
+  QRFSI12_Q%spknt(4) = 225.d0
+  QRFSI12_Q%spknt(5) = 300.d0
+  QRFSI12_Q%spknt(6) = 410.d0
+  QRFSI12_Q%spknt(7) = 530.d0
+  QRFSI12_Q%spknt(8) = 650.d0
 
 ! get the depths and 1/Q values of the reference model
   open(11,file=QRFSI12_ref,status='old',action='read',iostat=ier)
@@ -238,7 +238,7 @@
   else if(idoubling == IFLAG_OUTER_CORE_NORMAL) then
   !   print *,'QRFSI12: we are in the outer core'
      Qmu = 0.0d0
-  else !we are in the mantle
+  else ! we are in the mantle
     depth = 6371.-radius
 !   print *,'QRFSI12: we are in the mantle at depth',depth
     ifnd=0

Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/model_attenuation.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/model_attenuation.f90	2013-07-24 14:30:04 UTC (rev 22661)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/model_attenuation.f90	2013-07-24 15:33:12 UTC (rev 22662)
@@ -308,25 +308,25 @@
   ! uses "pure" 1D models including their 1D-crust profiles
   ! (uses USE_EXTERNAL_CRUSTAL_MODEL set to false)
   if(REFERENCE_1D_MODEL == REFERENCE_MODEL_PREM) then
-     AM_V%Qn = 12
+    AM_V%Qn = 12
   else if(REFERENCE_1D_MODEL == REFERENCE_MODEL_IASP91) then
-     AM_V%Qn = 12
+    AM_V%Qn = 12
   else if(REFERENCE_1D_MODEL == REFERENCE_MODEL_AK135F_NO_MUD) then
-     call define_model_ak135(.FALSE.,Mak135_V)
-     AM_V%Qn = NR_AK135F_NO_MUD
+    call define_model_ak135(.FALSE.,Mak135_V)
+    AM_V%Qn = NR_AK135F_NO_MUD
   else if(REFERENCE_1D_MODEL == REFERENCE_MODEL_1066A) then
-     call define_model_1066a(.FALSE.,M1066a_V)
-     AM_V%Qn = NR_1066A
+    call define_model_1066a(.FALSE.,M1066a_V)
+    AM_V%Qn = NR_1066A
   else if(REFERENCE_1D_MODEL == REFERENCE_MODEL_1DREF) then
-     call define_model_1dref(.FALSE.,Mref_V)
-     AM_V%Qn = NR_REF
+    call define_model_1dref(.FALSE.,Mref_V)
+    AM_V%Qn = NR_REF
   else if(REFERENCE_1D_MODEL == REFERENCE_MODEL_JP1D) then
-     AM_V%Qn = 12
+    AM_V%Qn = 12
   else if(REFERENCE_1D_MODEL == REFERENCE_MODEL_SEA1D) then
-     call define_model_sea1d(.FALSE.,SEA1DM_V)
-     AM_V%Qn = NR_SEA1D
+    call define_model_sea1d(.FALSE.,SEA1DM_V)
+    AM_V%Qn = NR_SEA1D
   else
-     call exit_MPI(myrank, 'Reference 1D Model Not recognized')
+    call exit_MPI(myrank, 'Reference 1D Model Not recognized')
   endif
 
   ! sets up attenuation storage (for all possible Qmu values defined in the 1D models)

Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/model_crust.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/model_crust.f90	2013-07-24 14:30:04 UTC (rev 22661)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/model_crust.f90	2013-07-24 15:33:12 UTC (rev 22662)
@@ -102,8 +102,8 @@
   type (model_crust_variables) CM_V
 ! model_crust_variables
 
-  double precision lat,lon,x,vp,vs,rho,moho
-  logical found_crust,elem_in_crust
+  double precision :: lat,lon,x,vp,vs,rho,moho
+  logical :: found_crust,elem_in_crust
 
   ! local parameters
   double precision :: h_sed,h_uc
@@ -173,19 +173,12 @@
 
   ! non-dimensionalize
   if (found_crust) then
-    scaleval = dsqrt(PI*GRAV*RHOAV)
-    vp = vp*1000.0d0/(R_EARTH*scaleval)
-    vs = vs*1000.0d0/(R_EARTH*scaleval)
-    rho = rho*1000.0d0/RHOAV
+    scaleval = ONE / ( R_EARTH_KM * dsqrt(PI*GRAV*RHOAV) )
+    vp = vp * scaleval
+    vs = vs * scaleval
+    rho = rho * 1000.0d0 / RHOAV
  endif
 
- ! checks moho value
- !moho = h_uc + thicks(6) + thicks(7)
- !if( moho /= thicks(NLAYERS_CRUST) ) then
- ! print*,'moho:',moho,thicks(NLAYERS_CRUST)
- ! print*,'  lat/lon/x:',lat,lon,x
- !endif
-
  ! No matter found_crust true or false, output moho thickness
  moho = (h_uc+thicks(6)+thicks(7))*1000.0d0/R_EARTH
 

Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/model_crustmaps.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/model_crustmaps.f90	2013-07-24 14:30:04 UTC (rev 22661)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/model_crustmaps.f90	2013-07-24 15:33:12 UTC (rev 22662)
@@ -80,8 +80,7 @@
   integer :: ier
 
   ! master reads in crust maps
-  if(myrank == 0) &
-    call read_general_crustmap(GC_V)
+  if(myrank == 0) call read_general_crustmap(GC_V)
 
   ! broadcasts values to all processes
   call MPI_BCAST(GC_V%thickness,180*360*CRUSTMAP_RESOLUTION*CRUSTMAP_RESOLUTION*NLAYERS_CRUSTMAP, &

Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/model_epcrust.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/model_epcrust.f90	2013-07-24 14:30:04 UTC (rev 22661)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/model_epcrust.f90	2013-07-24 15:33:12 UTC (rev 22662)
@@ -131,6 +131,7 @@
 
   subroutine model_epcrust(lat,lon,dep,vp,vs,rho,moho,found_crust,EPCRUST,elem_in_crust)
   implicit none
+
   include "constants.h"
 
   ! INPUT & OUTPUT

Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/model_heterogen_mantle.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/model_heterogen_mantle.f90	2013-07-24 14:30:04 UTC (rev 22661)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/model_heterogen_mantle.f90	2013-07-24 15:33:12 UTC (rev 22662)
@@ -99,7 +99,6 @@
   type (model_heterogen_m_variables) HMM
 ! model_heterogen_m_variables
 
-
   ! open heterogen.dat
   open(unit=10,file='./DATA/heterogen/heterogen.dat',access='direct',&
        form='formatted',recl=20,status='old',action='read',iostat=ier)

Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/model_jp3d.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/model_jp3d.f90	2013-07-24 14:30:04 UTC (rev 22661)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/model_jp3d.f90	2013-07-24 15:33:12 UTC (rev 22662)
@@ -439,9 +439,9 @@
 
   SUBROUTINE INPUT1(JP3DM_V)
 
-   implicit none
+  implicit none
 
-   include "constants.h"
+  include "constants.h"
 
 ! model_jp3d_variables
   type model_jp3d_variables
@@ -523,13 +523,20 @@
       RETURN
     END SUBROUTINE INPUT1
 
-      SUBROUTINE PUT1(NPX,NRX,NHX,PNX,RNX,HNX,VELXP)
-      integer :: NPX,NRX,NHX,K,I,J
-      double precision ::  VELXP(NPX,NRX,NHX), &
+!
+!-------------------------------------------------------------------------------------------------
+!
+  SUBROUTINE PUT1(NPX,NRX,NHX,PNX,RNX,HNX,VELXP)
+
+  implicit none
+
+  integer :: NPX,NRX,NHX,K,I,J
+  double precision ::  VELXP(NPX,NRX,NHX), &
                 PNX(NPX),RNX(NRX),HNX(NHX)
-      READ(2,110) (PNX(I),I=1,NPX)
-      READ(2,110) (RNX(I),I=1,NRX)
-      READ(2,120) (HNX(I),I=1,NHX)
+
+  READ(2,110) (PNX(I),I=1,NPX)
+  READ(2,110) (RNX(I),I=1,NRX)
+  READ(2,120) (HNX(I),I=1,NHX)
       DO K = 1,NHX
          DO I = 1,NPX
             READ(2,140) (VELXP(I,J,K),J=1,NRX)
@@ -544,6 +551,7 @@
 !---------------------------------------------------------------------------------------------
 !
   SUBROUTINE INPUT2(JP3DM_V)
+
   implicit none
 
   include "constants.h"
@@ -619,17 +627,17 @@
   type (model_jp3d_variables) JP3DM_V
 ! model_jp3d_variables
 
-      integer :: NP,NNR,I,J
-      READ(3,100)  NP,NNR
-      READ(3,110) (JP3DM_V%PN(I),I=1,NP)
-      READ(3,120) (JP3DM_V%RRN(I),I=1,NNR)
-      DO 1  I = NP,1,-1
+  integer :: NP,NNR,I,J
+  READ(3,100)  NP,NNR
+  READ(3,110) (JP3DM_V%PN(I),I=1,NP)
+  READ(3,120) (JP3DM_V%RRN(I),I=1,NNR)
+  DO 1  I = NP,1,-1
       READ(3,130) (JP3DM_V%DEPA(I,J),J=1,NNR)
 1     CONTINUE
-      DO 2  I = NP,1,-1
+  DO 2  I = NP,1,-1
       READ(3,130) (JP3DM_V%DEPB(I,J),J=1,NNR)
 2     CONTINUE
-      DO 3  I = NP,1,-1
+  DO 3  I = NP,1,-1
       READ(3,130) (JP3DM_V%DEPC(I,J),J=1,NNR)
 3     CONTINUE
 
@@ -649,6 +657,7 @@
   implicit none
 
   include "constants.h"
+
 ! model_jp3d_variables
   type model_jp3d_variables
     sequence

Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/model_s20rts.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/model_s20rts.f90	2013-07-24 14:30:04 UTC (rev 22661)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/model_s20rts.f90	2013-07-24 15:33:12 UTC (rev 22662)
@@ -34,7 +34,6 @@
 ! model, and that we use the PREM radial attenuation model when ATTENUATION is incorporated.
 !--------------------------------------------------------------------------------------------------
 
-
   subroutine model_s20rts_broadcast(myrank,S20RTS_V)
 
 ! standard routine to setup model
@@ -107,7 +106,7 @@
   ! local parameters
   integer :: k,l,m,ier
 
-  character(len=150) S20RTS, P12
+  character(len=150) :: S20RTS, P12
 
   call get_value_string(S20RTS, 'model.S20RTS', 'DATA/s20rts/S20RTS.dat')
   call get_value_string(P12, 'model.P12', 'DATA/s20rts/P12.dat')
@@ -259,10 +258,10 @@
   type (model_s20rts_variables) S20RTS_V
 ! model_s20rts_variables
 
+  ! local parameters
+  integer :: i,j
+  double precision :: qqwk(3,NK_20+1)
 
-  integer i,j
-  double precision qqwk(3,NK_20+1)
-
   S20RTS_V%spknt(1) = -1.00000d0
   S20RTS_V%spknt(2) = -0.78631d0
   S20RTS_V%spknt(3) = -0.59207d0



More information about the CIG-COMMITS mailing list