[4] | 1 | !> \file tab-litho-0.3.f90 |
---|
| 2 | !! Repartition des enfoncement en fonction de la rigidite de la lithosphere. |
---|
| 3 | !< |
---|
| 4 | |
---|
| 5 | !> SUBROUTINE: tab_litho |
---|
| 6 | !! \author ... |
---|
| 7 | !! \date ... |
---|
| 8 | !! @note Subroutine qui donne la repartition des enfoncements |
---|
| 9 | !! en fonction de la rigidite de la lithosphere. |
---|
| 10 | !! @note variables en entree |
---|
| 11 | !! - ROM masse volumique du manteau |
---|
| 12 | !! - RL rayon de rigidite relative |
---|
| 13 | !! - DL rigidite flexural |
---|
| 14 | !! |
---|
| 15 | !! @note variables en sortie |
---|
| 16 | !! - WE deflection due a une charge unitaire |
---|
| 17 | !! |
---|
| 18 | !! @note Used modules: |
---|
| 19 | !! @note - use module3D_phy |
---|
| 20 | !! @note - use param_phy_mod |
---|
| 21 | !! @note - use ISO_DECLAR |
---|
| 22 | !< |
---|
| 23 | |
---|
| 24 | subroutine tab_litho |
---|
| 25 | |
---|
| 26 | USE module3D_phy |
---|
| 27 | USE param_phy_mod |
---|
| 28 | USE ISO_DECLAR ! module de declaration pour l'isostasie |
---|
| 29 | |
---|
| 30 | |
---|
| 31 | implicit none |
---|
| 32 | |
---|
| 33 | |
---|
| 34 | INTEGER, PARAMETER :: nk=1000 |
---|
| 35 | REAL, dimension(0:nk) :: kei |
---|
| 36 | REAL :: stepk,AL,XL,DIST |
---|
| 37 | REAL :: som |
---|
| 38 | |
---|
| 39 | ! pour la lithosphere |
---|
| 40 | STEPK=100. |
---|
| 41 | AL=-RL*RL/(2.*PI*DL)*DX*DY |
---|
| 42 | |
---|
| 43 | ! fonction de kelvin |
---|
| 44 | ! lecture de la table kei qui est tous les 0.01 entre 0 et 10 |
---|
| 45 | ! STEPK=100=1/ecart |
---|
| 46 | |
---|
[10] | 47 | open(num_kelvin,file=trim(dirnameinp)//'../kelvin.res') |
---|
[4] | 48 | read(num_kelvin,*) |
---|
[53] | 49 | do K=0,NK |
---|
[4] | 50 | read(num_kelvin,*) XL,kei(K) |
---|
| 51 | end do |
---|
| 52 | close(num_kelvin) |
---|
| 53 | |
---|
| 54 | |
---|
| 55 | do I=-LBLOC,LBLOC |
---|
| 56 | do J=-LBLOC,LBLOC |
---|
| 57 | DIST=dx*sqrt(1.*(I*I+J*J)) |
---|
| 58 | XL=DIST/RL*STEPK |
---|
| 59 | K=int(XL) |
---|
| 60 | if ((K.gt.834).or.(DIST.gt.dx*LBLOC)) then |
---|
| 61 | We(I,J)=0. |
---|
| 62 | else |
---|
| 63 | We(I,J)=kei(k)+(kei(k+1)-kei(k))*(XL-K) |
---|
| 64 | if (K.eq.834) We(I,J)=max(We(I,J),0.) |
---|
| 65 | endif |
---|
| 66 | We(I,J)=We(I,J)*Al |
---|
| 67 | end do |
---|
| 68 | end do |
---|
| 69 | |
---|
| 70 | |
---|
| 71 | |
---|
| 72 | !open(55,file='iso-cat') |
---|
| 73 | !write(55,*) 'WE avec LBLOC sur 480 km' |
---|
| 74 | !do i=-lbloc,lbloc |
---|
| 75 | ! write(55,'25(i4,1x)') (int(WE(i,j)*1.e08),j=-lbloc,lbloc) |
---|
| 76 | !end do |
---|
| 77 | !close(55) |
---|
| 78 | |
---|
| 79 | ! sommation des elements de WE dans som pour normaliser |
---|
| 80 | ! la somme des enfoncements a 1 |
---|
| 81 | som = SUM(WE(:,:)) |
---|
| 82 | |
---|
| 83 | |
---|
| 84 | ! normalisation |
---|
| 85 | WE(:,:) = WE(:,:)/(som*ROM*g) |
---|
| 86 | |
---|
| 87 | |
---|
| 88 | |
---|
| 89 | ! ATTENTION DESSOUS VALABLE SEULEMENT EN 40 km |
---|
| 90 | |
---|
| 91 | ! la partie ci-dessous doit etre mise dans une routine specifique Philippe |
---|
| 92 | ! pour faire exactement la meme isostasie que Philippe |
---|
| 93 | !if ((icouple.eq.2).or.(icouple.eq.4)) then |
---|
| 94 | |
---|
| 95 | ! LBLOC=10 |
---|
| 96 | |
---|
| 97 | ! do I=-LBLOC,LBLOC |
---|
| 98 | ! do J=-LBLOC,LBLOC |
---|
| 99 | ! We(i,j)=0. ! mise a 0 partout |
---|
| 100 | ! end do |
---|
| 101 | ! end do |
---|
| 102 | ! open(num_phil,file='../INPUT/iso-philippe') |
---|
| 103 | |
---|
| 104 | ! read(num_phil,*) |
---|
| 105 | ! read(num_phil,*) |
---|
| 106 | ! read(num_phil,*) |
---|
| 107 | ! read(num_phil,*) |
---|
| 108 | ! do J=-10,10 |
---|
| 109 | ! read(num_phil,*)(We(i,j),i=-10,10) |
---|
| 110 | ! end do |
---|
| 111 | ! close(num_phil) |
---|
| 112 | ! do i=-LBLOC,LBLOC |
---|
| 113 | ! do j=-LBLOC,LBLOC |
---|
| 114 | ! we(i,j)=we(i,j)/11030/g/rom |
---|
| 115 | ! end do |
---|
| 116 | ! end do |
---|
| 117 | !endif |
---|
| 118 | |
---|
| 119 | return |
---|
| 120 | end subroutine tab_litho |
---|