[cig-commits] [commit] devel: Modif HEX27 (7d2d571)
cig_noreply at geodynamics.org
cig_noreply at geodynamics.org
Thu Nov 20 10:08:12 PST 2014
Repository : https://github.com/geodynamics/specfem3d
On branch : devel
Link : https://github.com/geodynamics/specfem3d/compare/b43c63c4cea034e3e6f647cebc591fb4c2ab3894...fb5ff57e01f82037a9b22ebd3c3c184e9e4f14f6
>---------------------------------------------------------------
commit 7d2d5716e892a9483a247e2157e6dba4a4888ed8
Author: Clément Durochat <c.durochat at gmail.com>
Date: Wed Nov 12 11:09:30 2014 +0100
Modif HEX27
>---------------------------------------------------------------
7d2d5716e892a9483a247e2157e6dba4a4888ed8
src/meshfem3D/earth_chunk_HEX27_Mesher.f90 | 595 +----------------------------
src/meshfem3D/earth_chunk_ReadIasp91.f90 | 5 +-
2 files changed, 4 insertions(+), 596 deletions(-)
diff --git a/src/meshfem3D/earth_chunk_HEX27_Mesher.f90 b/src/meshfem3D/earth_chunk_HEX27_Mesher.f90
index 2638757..4b276be 100644
--- a/src/meshfem3D/earth_chunk_HEX27_Mesher.f90
+++ b/src/meshfem3D/earth_chunk_HEX27_Mesher.f90
@@ -809,597 +809,4 @@
end subroutine earth_chunk_HEX27_Mesher
-!=======================================================================================================
-!
-!=======================================================================================================
-
-! compute the Euler angles and the associated rotation matrix
-
- subroutine euler_angles(rotation_matrix,CENTER_LONGITUDE_IN_DEGREES,CENTER_LATITUDE_IN_DEGREES,GAMMA_ROTATION_AZIMUTH)
-
- implicit none
-
- !include "constants.h"
-
- double precision rotation_matrix(3,3)
- double precision CENTER_LONGITUDE_IN_DEGREES,CENTER_LATITUDE_IN_DEGREES,GAMMA_ROTATION_AZIMUTH
-
- double precision alpha,beta,gamma
- double precision sina,cosa,sinb,cosb,sing,cosg
-
- double precision DEGREES_TO_RADIANS
-
- DEGREES_TO_RADIANS = 3.141592653589793d0/180.d0
-
-
-! compute colatitude and longitude and convert to radians
- alpha = CENTER_LONGITUDE_IN_DEGREES * DEGREES_TO_RADIANS
- beta = (90.0d0 - CENTER_LATITUDE_IN_DEGREES) * DEGREES_TO_RADIANS
- gamma = GAMMA_ROTATION_AZIMUTH * DEGREES_TO_RADIANS
-
- sina = dsin(alpha)
- cosa = dcos(alpha)
- sinb = dsin(beta)
- cosb = dcos(beta)
- sing = dsin(gamma)
- cosg = dcos(gamma)
-
-! define rotation matrix
- rotation_matrix(1,1) = cosg*cosb*cosa-sing*sina
- rotation_matrix(1,2) = -sing*cosb*cosa-cosg*sina
- rotation_matrix(1,3) = sinb*cosa
- rotation_matrix(2,1) = cosg*cosb*sina+sing*cosa
- rotation_matrix(2,2) = -sing*cosb*sina+cosg*cosa
- rotation_matrix(2,3) = sinb*sina
- rotation_matrix(3,1) = -cosg*sinb
- rotation_matrix(3,2) = sing*sinb
- rotation_matrix(3,3) = cosb
-
- end subroutine euler_angles
-
-!=======================================================================================================
-!
-!=======================================================================================================
-
- subroutine write_gllz_points(xstore,ystore,zstore,NGLLX,NGLLY,NGLLZ,current_layer,nel_depth,ilayer,iz,Ndepth)
-
- implicit none
-
- integer NGLLX,NGLLY,NGLLZ,nel_depth,iz,Ndepth
- double precision xstore(NGLLX,NGLLY,NGLLZ),ystore(NGLLX,NGLLY,NGLLZ),zstore(NGLLX,NGLLY,NGLLZ)
- double precision profondeur
- integer current_layer(0:nel_depth-1),ilayer,k
- !write(*,*) ilayer, current_layer(iz)
- !profondeur = dsqrt(xstore(1,1,k)**2 + ystore(1,1,k)**2 + (zstore(1,1,k) )**2 )
- !write(27,*) profondeur/1000., ilayer
- if (ilayer == current_layer(iz)) then
- do k=2,NGLLZ
- profondeur = dsqrt(xstore(1,1,k)**2 + ystore(1,1,k)**2 + (zstore(1,1,k) )**2 )
- write(27,*) profondeur/1000., ilayer-1,1
- Ndepth=Ndepth+1
- enddo
- else ! new layer
-
- k=1
- profondeur = dsqrt(xstore(1,1,k)**2 + ystore(1,1,k)**2 + (zstore(1,1,k) )**2 )
- if (ilayer==0) then
- ilayer = current_layer(iz)
- write(27,*) profondeur/1000., ilayer-1,1
- Ndepth=Ndepth+1
- else
- ilayer = current_layer(iz)
- write(27,*) profondeur/1000., ilayer-1,-1
- Ndepth=Ndepth+1
- endif
- do k=2,NGLLZ ! on duplique le dernier point
- profondeur = dsqrt(xstore(1,1,k)**2 + ystore(1,1,k)**2 + (zstore(1,1,k) )**2 )
- write(27,*) profondeur/1000., ilayer-1,1
- Ndepth=Ndepth+1
- enddo
-
-
- endif
-
- end subroutine write_gllz_points
-
-!=======================================================================================================
-!
-!=======================================================================================================
-
- subroutine write_recdepth_dsm(Ndepth,R_EARTH,MESH)
-
- implicit none
-
- integer Ndepth,i
- double precision R_EARTH,prof
- double precision, allocatable :: z(:)
- integer, allocatable :: zindex(:),ziflag(:)
- integer ilayer,flag
- character(len=10) MESH
-
-
-
- open(27,file=trim(MESH)//'.recdepth')
-
- allocate(zindex(Ndepth),ziflag(Ndepth))
- allocate(z(Ndepth))
-
- do i=1,Ndepth
- read(27,*) prof,ilayer,flag
- z(Ndepth-i+1)=R_EARTH/1000.d0-prof
- zindex(Ndepth-i+1)=ilayer
- ziflag(Ndepth-i+1)=flag
- enddo
-
- close(27)
- open(27,file=trim(MESH)//'recdepth')
-
- write(27,*) Ndepth
- i=1
- write(27,*) z(i),zindex(i),ziflag(i)
- do i=2,Ndepth-1
- if (ziflag(i-1) == -1 ) then
- write(27,*) z(i),zindex(i),-1
- else
- write(27,*) z(i),zindex(i),1
- endif
- enddo
- i=Ndepth
- write(27,*) z(i),zindex(i),ziflag(i)
-
-end subroutine write_recdepth_dsm
-
-!=======================================================================================================
-!
-!=======================================================================================================
-
- subroutine write_stxmin(xstore,ystore,zstore,NGLLX,NGLLY,NGLLZ,rotation_matrix,test)
-
- implicit none
-
- integer NDIM,NGLLX,NGLLY,NGLLZ,jgll,i,j,NGLLY_eff
- double precision xstore(NGLLX,NGLLY,NGLLZ),ystore(NGLLX,NGLLY,NGLLZ),zstore(NGLLX,NGLLY,NGLLZ)
- double precision rotation_matrix(3,3)
- double precision vector_ori(3),vector_rotated(3)
- double precision rayon,x,y,z,deg2rad,long,lati
- logical test
-
- deg2rad=3.141592653589793d0/180.d0
- NDIM=3
-
- if (test) then
- NGLLY_eff = NGLLY
- else
- NGLLY_eff = NGLLY - 1
- endif
-
- do jgll=1,NGLLY_eff
- vector_ori(1)=xstore(1,jgll,NGLLZ)
- vector_ori(2)=ystore(1,jgll,NGLLZ)
- vector_ori(3)=zstore(1,jgll,NGLLZ)
-
- do i = 1,NDIM
- vector_rotated(i) = 0.d0
- do j = 1,NDIM
- vector_rotated(i) = vector_rotated(i) + rotation_matrix(i,j)*vector_ori(j)
- enddo
- enddo
- x=vector_rotated(1);y=vector_rotated(2);z=vector_rotated(3)
- rayon = dsqrt(vector_rotated(1)**2 + vector_rotated(2)**2 + vector_rotated(3)**2)
-
- long=atan2(y,x)
- lati=asin(z/rayon)
-
- ! passage de geocentique à géographique
- !!theta = PI/2.D0 - lati
- ! convert the geocentric colatitude to a geographic colatitude
- !!colat = PI/2.0d0 - datan(1.006760466d0*dcos(theta)/dmax1(TINYVAL,dsin(theta)))
- !!lati = PI/2.0d0 - colat
-
- !write(28,*) xstore(1,jgll,NGLLZ), ystore(1,jgll,NGLLZ), zstore(1,jgll,NGLLZ)!x,y !long/deg2rad,lati/deg2rad
- write(28,*) long/deg2rad,lati/deg2rad !,rayon/1000
- !write(38,'()') 1,(NGLLY-1)*jy_elm+jgll
- write(49,*)
- write(49,*) vector_ori(:)
- write(49,*) vector_rotated(:)
-
- enddo
-
- end subroutine write_stxmin
-
-!=======================================================================================================
-!
-!=======================================================================================================
-
- subroutine write_stxmax(xstore,ystore,zstore,NGLLX,NGLLY,NGLLZ,rotation_matrix,test)
-
- implicit none
-
- integer NDIM,NGLLX,NGLLY,NGLLZ,jgll,i,j,NGLLY_eff
- double precision xstore(NGLLX,NGLLY,NGLLZ),ystore(NGLLX,NGLLY,NGLLZ),zstore(NGLLX,NGLLY,NGLLZ)
- double precision rotation_matrix(3,3)
- double precision vector_ori(3),vector_rotated(3)
- double precision rayon,x,y,z,deg2rad,long,lati
- logical test
-
- if (test) then
- NGLLY_eff = NGLLY
- else
- NGLLY_eff = NGLLY - 1
- endif
-
- deg2rad=3.141592653589793d0/180.d0
- NDIM=3
-
- do jgll=1,NGLLY_eff
- vector_ori(1)=xstore(NGLLX,jgll,NGLLZ)
- vector_ori(2)=ystore(NGLLX,jgll,NGLLZ)
- vector_ori(3) =zstore(NGLLX,jgll,NGLLZ)
-
- do i = 1,NDIM
- vector_rotated(i) = 0.d0
- do j = 1,NDIM
- vector_rotated(i) = vector_rotated(i) + rotation_matrix(i,j)*vector_ori(j)
- enddo
- enddo
- x=vector_rotated(1);y=vector_rotated(2);z=vector_rotated(3)
- rayon = dsqrt(vector_rotated(1)**2 + vector_rotated(2)**2 + vector_rotated(3)**2)
-
- long=atan2(y,x)
- lati=asin(z/rayon)
-
- ! passage de geocentique à géographique
- !!theta = PI/2.D0 - lati
- ! convert the geocentric colatitude to a geographic colatitude
- !!colat = PI/2.0d0 - datan(1.006760466d0*dcos(theta)/dmax1(TINYVAL,dsin(theta)))
- !!lati = PI/2.0d0 - colat
-
- !write(28,*) xstore(1,jgll,NGLLZ), ystore(1,jgll,NGLLZ), zstore(1,jgll,NGLLZ)!x,y !long/deg2rad,lati/deg2rad
- write(29,*) long/deg2rad,lati/deg2rad !,rayon/1000
- enddo
-
- end subroutine write_stxmax
-
-!=======================================================================================================
-!
-!=======================================================================================================
-
- subroutine write_stymin(xstore,ystore,zstore,NGLLX,NGLLY,NGLLZ,rotation_matrix,test)
-
- implicit none
-
- integer NDIM,NGLLX,NGLLY,NGLLZ,jgll,i,j,NGLLX_eff
- double precision xstore(NGLLX,NGLLY,NGLLZ),ystore(NGLLX,NGLLY,NGLLZ),zstore(NGLLX,NGLLY,NGLLZ)
- double precision rotation_matrix(3,3)
- double precision vector_ori(3),vector_rotated(3)
- double precision rayon,x,y,z,deg2rad,long,lati
- logical test
-
- deg2rad=3.141592653589793d0/180.d0
- NDIM=3
-
- if (test) then
- NGLLX_eff = NGLLX
- else
- NGLLX_eff = NGLLX - 1
- endif
-
- do jgll=1,NGLLX_eff
- vector_ori(1)=xstore(jgll,1,NGLLZ)
- vector_ori(2)=ystore(jgll,1,NGLLZ)
- vector_ori(3) =zstore(jgll,1,NGLLZ)
-
- do i = 1,NDIM
- vector_rotated(i) = 0.d0
- do j = 1,NDIM
- vector_rotated(i) = vector_rotated(i) + rotation_matrix(i,j)*vector_ori(j)
- enddo
- enddo
- x=vector_rotated(1);y=vector_rotated(2);z=vector_rotated(3)
- rayon = dsqrt(vector_rotated(1)**2 + vector_rotated(2)**2 + vector_rotated(3)**2)
-
- long=atan2(y,x)
- lati=asin(z/rayon)
-
- ! passage de geocentique à géographique
- !!theta = PI/2.D0 - lati
- ! convert the geocentric colatitude to a geographic colatitude
- !!colat = PI/2.0d0 - datan(1.006760466d0*dcos(theta)/dmax1(TINYVAL,dsin(theta)))
- !!lati = PI/2.0d0 - colat
-
- !write(28,*) xstore(1,jgll,NGLLZ), ystore(1,jgll,NGLLZ), zstore(1,jgll,NGLLZ)!x,y !long/deg2rad,lati/deg2rad
- write(30,*) long/deg2rad,lati/deg2rad !,rayon/1000
- enddo
-
- end subroutine write_stymin
-
-!=======================================================================================================
-!
-!=======================================================================================================
-
- subroutine write_stymax(xstore,ystore,zstore,NGLLX,NGLLY,NGLLZ,rotation_matrix,test)
-
- implicit none
-
- integer NDIM,NGLLX,NGLLY,NGLLZ,jgll,i,j,NGLLX_eff
- double precision xstore(NGLLX,NGLLY,NGLLZ),ystore(NGLLX,NGLLY,NGLLZ),zstore(NGLLX,NGLLY,NGLLZ)
- double precision rotation_matrix(3,3)
- double precision vector_ori(3),vector_rotated(3)
- double precision rayon,x,y,z,deg2rad,long,lati
- logical test
-
- if (test) then
- NGLLX_eff = NGLLX
- else
- NGLLX_eff = NGLLX - 1
- endif
-
- deg2rad=3.141592653589793d0/180.d0
- NDIM=3
-
- do jgll=1,NGLLX_eff
- vector_ori(1)=xstore(jgll,NGLLY,NGLLZ)
- vector_ori(2)=ystore(jgll,NGLLY,NGLLZ)
- vector_ori(3) =zstore(jgll,NGLLY,NGLLZ)
-
- do i = 1,NDIM
- vector_rotated(i) = 0.d0
- do j = 1,NDIM
- vector_rotated(i) = vector_rotated(i) + rotation_matrix(i,j)*vector_ori(j)
- enddo
- enddo
- x=vector_rotated(1);y=vector_rotated(2);z=vector_rotated(3)
- rayon = dsqrt(vector_rotated(1)**2 + vector_rotated(2)**2 + vector_rotated(3)**2)
-
- long=atan2(y,x)
- lati=asin(z/rayon)
-
- ! passage de geocentique à géographique
- !!theta = PI/2.D0 - lati
- ! convert the geocentric colatitude to a geographic colatitude
- !!colat = PI/2.0d0 - datan(1.006760466d0*dcos(theta)/dmax1(TINYVAL,dsin(theta)))
- !!lati = PI/2.0d0 - colat
-
- !write(28,*) xstore(1,jgll,NGLLZ), ystore(1,jgll,NGLLZ), zstore(1,jgll,NGLLZ)!x,y !long/deg2rad,lati/deg2rad
- write(31,*) long/deg2rad,lati/deg2rad !,rayon/1000
- enddo
-
- end subroutine write_stymax
-
-!=======================================================================================================
-!
-!=======================================================================================================
-
- subroutine store_zmin_points(xstore,ystore,zstore,NGLLX,NGLLY,NGLLZ,rotation_matrix,&
- lon_zmin,lat_zmin,nlon_dsm,nlat_dsm,ilon,ilat)
-
- implicit none
-
- integer NDIM,NGLLX,NGLLY,NGLLZ,igll,jgll,i,j
- integer ilon,ilat,iglob,jglob,nlat_dsm,nlon_dsm
- double precision xstore(NGLLX,NGLLY,NGLLZ),ystore(NGLLX,NGLLY,NGLLZ),zstore(NGLLX,NGLLY,NGLLZ)
- double precision rotation_matrix(3,3)
- double precision vector_ori(3),vector_rotated(3)
- double precision rayon,x,y,z,deg2rad,long,lati
- double precision lon_zmin(nlon_dsm,nlat_dsm),lat_zmin(nlon_dsm,nlat_dsm)
-
-
- deg2rad=3.141592653589793d0/180.d0
- NDIM=3
-
- do jgll=1,NGLLY
- do igll=1,NGLLX
- vector_ori(1)=xstore(igll,jgll,1)
- vector_ori(2)=ystore(igll,jgll,1)
- vector_ori(3) =zstore(igll,jgll,1)
-
- do i = 1,NDIM
- vector_rotated(i) = 0.d0
- do j = 1,NDIM
- vector_rotated(i) = vector_rotated(i) + rotation_matrix(i,j)*vector_ori(j)
- enddo
- enddo
- x=vector_rotated(1);y=vector_rotated(2);z=vector_rotated(3)
- rayon = dsqrt(vector_rotated(1)**2 + vector_rotated(2)**2 + vector_rotated(3)**2)
-
- long=atan2(y,x)
- lati=asin(z/rayon)
-
- ! passage de geocentique à géographique
- !!theta = PI/2.D0 - lati
- ! convert the geocentric colatitude to a geographic colatitude
- !!colat = PI/2.0d0 - datan(1.006760466d0*dcos(theta)/dmax1(TINYVAL,dsin(theta)))
- !!lati = PI/2.0d0 - colat
-
- !write(28,*) xstore(1,jgll,NGLLZ), ystore(1,jgll,NGLLZ), zstore(1,jgll,NGLLZ)!x,y !long/deg2rad,lati/deg2rad
- !write(31,*) long/deg2rad,lati/deg2rad !,rayon/1000
- iglob=(ilon)*(NGLLX-1)+igll
- jglob=(ilat)*(NGLLY-1)+jgll
- lon_zmin(iglob,jglob)= long/deg2rad
- lat_zmin(iglob,jglob)= lati/deg2rad
- !write(32,'(3f20.10)') xstore(igll,jgll,1)/1000.d0, ystore(igll,jgll,1)/1000.d0,zstore(igll,jgll,1)/1000.d0
- !write(32,*) xstore(igll,jgll,NGLLZ), ystore(igll,igll,NGLLZ),zstore(igll,jgll,NGLLZ)
- enddo
- enddo
-
- end subroutine store_zmin_points
-
-!=======================================================================================================
-!
-!=======================================================================================================
-
- subroutine write_stzmin(x,y,nx,ny,MESH)
-
- implicit none
-
- integer i,j,nx,ny
- double precision x(nx,ny),y(nx,ny)
- character(len=10) MESH
-
- open(27,file=trim(MESH)//'stzmin')
- write(27,*) nx*ny
- do j=1,ny
- do i=1,nx
- write(27,*) x(i,j),y(i,j)
- enddo
- enddo
- close(27)
-
- end subroutine write_stzmin
-
-!=======================================================================================================
-!
-!=======================================================================================================
-
- subroutine write_Igm_file(iunit,ispec2D,NGLL1,NGLL2,ie,je,js,il)
-
- implicit none
-
- integer iunit,ispec2D,NGLL1,NGLL2,ie,je,js,il
- integer i,j
- do j=1,NGLL2
- do i=1,NGLL1
- write(iunit,*) i,j,ispec2D,(NGLL1-1)*ie+i,(NGLL2-1)*je+j+js,il
- enddo
- enddo
-
- end subroutine write_Igm_file
-
-!=======================================================================================================
-!
-!=======================================================================================================
-
- subroutine compute_rotation_matrix(rotation_matrix, lon_center_chunk,lat_center_chunk, chunk_azi)
-
- implicit none
-
- double precision rotation_matrix(3,3),lon_center_chunk,lat_center_chunk, chunk_azi
- double precision R0(3,3),R1(3,3),R2(3,3),axe_rotation(3),R00(3,3)
-
- ! je met le chunk en 0,0
- axe_rotation(1)=0.d0; axe_rotation(2)=1.d0; axe_rotation(3)=0.d0
- call rotation_matrix_axe(R00,axe_rotation,90.d0) ! je ramene le chunk en (0,0)
- ! rotation de l'azimuth du chunk
- axe_rotation(1)=1.d0; axe_rotation(2)=0.d0; axe_rotation(3)=0.d0
- call rotation_matrix_axe(R0,axe_rotation,90.-chunk_azi)
- ! on met le chunk a la bonne latitude
- axe_rotation(1)=0.d0; axe_rotation(2)=-1.d0; axe_rotation(3)=0.d0
- call rotation_matrix_axe(R1,axe_rotation,lat_center_chunk)
- ! on met le chunk a la bonne longitude
- axe_rotation(1)=0.d0; axe_rotation(2)=0.d0; axe_rotation(3)=1.d0
- call rotation_matrix_axe(R2,axe_rotation, lon_center_chunk)
- ! rotation resultante
- call compose4matrix(rotation_matrix,R00,R0,R1,R2)
-
- end subroutine compute_rotation_matrix
-
-!=======================================================================================================
-!
-! ROUTINES POUR FAIRE DES ROTATIONS 3D ET DIVERS CHANGEMENTS DE REPERES
-!
-! Vadim Monteiller Mars 2013
-!
-!-------------------------------------------------------------------------------
-! matrice de rotation 3D d'axe "axe" et d'angle theta (d°)
-! cette matrice est en complexe
-!
-!=======================================================================================================
-!
- subroutine rotation_matrix_axe(R,axe,theta)
-
- implicit none
-
- double precision axe(3),theta,pi,deg2rad
- double precision R(3,3)
- double precision c,s,ux,uy,uz,norme_axe
-
- pi=3.1415926535897932d0
- deg2rad = pi / 180.d0
- ! on normalise l'axe
- norme_axe=dsqrt(axe(1)**2 + axe(2)**2 + axe(3)**2)
-
- ! composantes de l'axe
- ux=axe(1)/norme_axe
- uy=axe(2)/norme_axe
- uz=axe(3)/norme_axe
-
- ! on calcule le cos et sin
- c=dcos(deg2rad * theta);s=dsin(deg2rad * theta)
-
- ! matrice de rotation complexe
- R(1,1)=(ux**2 + (1.d0-ux**2)*c)
- R(1,2)=(ux*uy*(1.d0-c)-uz*s)
- R(1,3)=(ux*uy*(1.d0-c)+uy*s)
-
- R(2,1)=(ux*uy*(1.d0-c)+uz*s)
- R(2,2)=(uy**2+(1.d0-uy**2)*c)
- R(2,3)=(uy*uz*(1.d0-c)-ux*s)
-
- R(3,1)=(ux*uz*(1.d0-c)-uy*s)
- R(3,2)=(uy*uz*(1.d0-c)+ux*s)
- R(3,3)=(uz**2+(1.d0-uz**2)*c)
-
- write(49,*) ' MATRICE ROTATION '
- write(49,*) R(1,:)
- write(49,*) R(2,:)
- write(49,*) R(3,:)
- write(49,*)
-
- end subroutine rotation_matrix_axe
-
-!=======================================================================================================
-!
-! R=R2*R1*R0
-!
-!=======================================================================================================
-
- subroutine compose4matrix(R,R00,R0,R1,R2)
-
- implicit none
-
- double precision R(3,3),R0(3,3),R1(3,3),R2(3,3),R00(3,3),Rtmp(3,3)
- integer i,j,k
-
-
- R(:,:)=0.d0
- ! multiplication R=R0*R00
- do j=1,3
- do i=1,3
- do k=1,3
- R(i,j)=R(i,j) + R0(i,k)*R00(k,j)
- enddo
- enddo
- enddo
-
- ! multiplication R=R1*R
- Rtmp=R
- R(:,:)=0.d0
- do j=1,3
- do i=1,3
- do k=1,3
- R(i,j)=R(i,j) + R1(i,k)*Rtmp(k,j)
- enddo
- enddo
- enddo
-
- ! multiplication R=R2*R
- Rtmp=R
- R(:,:)=0.d0
- do j=1,3
- do i=1,3
- do k=1,3
- R(i,j)=R(i,j) + R2(i,k)*Rtmp(k,j)
- enddo
- enddo
- enddo
-
- write(49,*) ' MATRICE ROTATION COMPLETE '
- write(49,*) R(1,:)
- write(49,*) R(2,:)
- write(49,*) R(3,:)
- write(49,*)
-
- end subroutine compose4matrix
-
-!------------------------------------------------------------------------------
-! rotation pour passer d'un repere local a un autre
+!=======================================================================================================!
diff --git a/src/meshfem3D/earth_chunk_ReadIasp91.f90 b/src/meshfem3D/earth_chunk_ReadIasp91.f90
index 7762da7..20c0994 100644
--- a/src/meshfem3D/earth_chunk_ReadIasp91.f90
+++ b/src/meshfem3D/earth_chunk_ReadIasp91.f90
@@ -141,8 +141,9 @@
end subroutine Read_dsm_model
-!
-!===========================================================================
+!===========================================================================!
+
+!===========================================================================!
!
subroutine Lyfnd(r,rb,n,i)
More information about the CIG-COMMITS
mailing list