1 | !> \file litho-0.4.f90 |
---|
2 | !! Donne la repartition des enfoncements |
---|
3 | !! en fonction de la rigidite de la lithosphere. |
---|
4 | !< |
---|
5 | |
---|
6 | !> SUBROUTINE: litho() |
---|
7 | !!Petit routine qui donne la repartition des enfoncements |
---|
8 | !!en fonction de la rigidite de la lithosphere. |
---|
9 | !! \author .... |
---|
10 | !! \date 10 septembre 2009 |
---|
11 | !! @note Cette version utilise la syntaxe f90 pour calculer W1 |
---|
12 | !! - avantage : plus facile a lire |
---|
13 | !! pas de pb d'execution avac ifor 10 |
---|
14 | !! - inconvenient : pourrait perdre la symetrie. |
---|
15 | !! @note - En entree |
---|
16 | !! @note WE(-LBLOC:LBLOC,-LBLOC:LBLOC) : deflection due a une charge unitaire |
---|
17 | !! defini dans tab-litho |
---|
18 | !! @note LBLOC : relie a la distance : distance en noeud autour de laquelle |
---|
19 | !! la flexure de la lithosphere est calculee |
---|
20 | !! @note CHARGE(1-LBLOC:NX+LBLOC,1-LBLOC:NY+LBLOC) : poids par unite de surface |
---|
21 | !! (unite ?) au temps time, calcule avant 'appel a litho |
---|
22 | !! dans taubed ou initial2 |
---|
23 | !! @note - En sortie |
---|
24 | !! @note W1(nx,ny) : enfoncement courant |
---|
25 | !! @note - tableaux de travail. |
---|
26 | !! @note WLOC(-LBLOC:LBLOC,-LBLOC:LBLOC) : enfoncement au point i,j du a la charge |
---|
27 | !! au point IP,JP |
---|
28 | !! @note LPX et LPY sont les coordonnees entre (-LBLOC:LBLOC,-LBLOC:LBLOC) autour |
---|
29 | !! du point I,J pour calculer tous les WLOC |
---|
30 | !! @note IP et JP sont les coordonnees de CHARGE dans un repere : |
---|
31 | !! (1-LBLOC:NX+LBLOC,1-LBLOC:NY+LBLOC) |
---|
32 | !! @note Used modules: |
---|
33 | !! @note - use module3D_phy |
---|
34 | !! @note - use ISO_DECLAR |
---|
35 | !! @note - use param_phy_mod |
---|
36 | !< |
---|
37 | |
---|
38 | |
---|
39 | subroutine litho |
---|
40 | !$ USE OMP_LIB |
---|
41 | USE module3D_phy |
---|
42 | USE param_phy_mod |
---|
43 | USE ISO_DECLAR ! module de declaration des variables specifiques a l'isostasie |
---|
44 | |
---|
45 | implicit none |
---|
46 | |
---|
47 | INTEGER :: SOM1,SOM2 |
---|
48 | REAL, dimension(:,:), allocatable :: WLOC |
---|
49 | |
---|
50 | |
---|
51 | !----- allocation de WLOC et de croix ----------- |
---|
52 | |
---|
53 | if (.not.allocated(WLOC)) THEN |
---|
54 | allocate(WLOC(-LBLOC:LBLOC,-LBLOC:LBLOC),stat=err) |
---|
55 | if (err/=0) then |
---|
56 | print *,"Erreur a l'allocation du tableau WLOC",err |
---|
57 | stop 4 |
---|
58 | end if |
---|
59 | end if |
---|
60 | |
---|
61 | !----- fin de l'allocation -------------- |
---|
62 | |
---|
63 | ! calcul de la deflexion par sommation des voisins appartenant |
---|
64 | ! au bloc de taille 2*LBLOC |
---|
65 | som1=0. |
---|
66 | som2=0. |
---|
67 | |
---|
68 | ! On somme aussi les contributions des points exterieurs au domaine |
---|
69 | ! lorsque la charge est due a l'ocean. On suppose alors que |
---|
70 | ! ces points ont la meme charge que les limites |
---|
71 | !$OMP PARALLEL |
---|
72 | !$OMP DO PRIVATE(Wloc) |
---|
73 | boucleij : do J=1,NY |
---|
74 | do I=1,NX |
---|
75 | |
---|
76 | |
---|
77 | W1(I,J)=0. |
---|
78 | ! ii=0 |
---|
79 | |
---|
80 | ! Wloc : enfoncement centre sur ij du a chaque charges autour |
---|
81 | Wloc(:,:) = We(:,:) * charge(i-lbloc:i+lbloc,j-lbloc:j+lbloc) |
---|
82 | |
---|
83 | |
---|
84 | W1(i,j) = sum (Wloc(:,:)) ! effet en ij de toutes les charges autour (distance < lbloc) |
---|
85 | |
---|
86 | ! ------------------------ FIN DE L'INTEGRATION SUR LE PAVE LBLOC |
---|
87 | ! som1=som1+W1(I,J) |
---|
88 | ! som2=som2-charge(I,J)/ROMG |
---|
89 | |
---|
90 | end do |
---|
91 | end do boucleij |
---|
92 | !$OMP END DO |
---|
93 | !$OMP END PARALLEL |
---|
94 | |
---|
95 | ! write(6,*)'enfoncement total avec rigidite de la lithosphere:', |
---|
96 | ! & som1 |
---|
97 | ! write(6,*)'enfoncement total avec isostasie locale :', |
---|
98 | ! & som2 |
---|
99 | |
---|
100 | |
---|
101 | end subroutine litho |
---|
102 | |
---|