|
Oasis3 4.0.2
|
|
Oasis3 4.0.2
|
00001 SUBROUTINE conserv (pfldb, ksizb, kmskb, psurfb, 00002 $ pflda, ksiza, kmska, psurfa, 00003 $ psgrb, psgra, cdmet) 00004 C**** 00005 C ***************************** 00006 C * OASIS ROUTINE - LEVEL 3 * 00007 C * ------------- ------- * 00008 C ***************************** 00009 C 00010 C**** *conserv* - Flux conservation routine 00011 C 00012 C Purpose: 00013 C ------- 00014 C Use global lagrange multiplier to insure conservation 00015 C 00016 C** Interface: 00017 C --------- 00018 C *CALL* *conserv (pfldb, ksizb, kmskb, psurfb, 00019 C pflda, ksiza, kmska, psurfa, cdmet)* 00020 C 00021 C Input: 00022 C ----- 00023 C pflda : field on target grid (real 1D) 00024 C pfldb : field on source grid (real 1D) 00025 C kmska : mask for target grid (integer 1D) 00026 C kmskb : mask for source grid (integer 1D) 00027 C psurfa : surfaces for target grid meshes (real 1D) 00028 C psurfb : surfaces for source grid meshes (real 1D) 00029 C ksizb : source arrays size (integer) 00030 C ksiza : target arrays size (integer) 00031 C psgrb : work array (real 1D) 00032 C psgra : work array (real 1D) 00033 C cdmet : conservation method (character string) 00034 C 00035 C Output: 00036 C ------ 00037 C pflda : field on target grid with conservation (real 1D) 00038 C 00039 C Workspace: 00040 C --------- 00041 C None 00042 C 00043 C Externals: 00044 C --------- 00045 C None 00046 C 00047 C Reference: 00048 C --------- 00049 C See OASIS manual (2000) 00050 C 00051 C History: 00052 C ------- 00053 C Version Programmer Date Description 00054 C ------- ---------- ---- ----------- 00055 C 1.0 L. Terray 94/01/01 created 00056 C 2.0 L. Terray 95/10/01 modified: new structure 00057 C 2.1 L. Terray 96/09/01 modified: printing 00058 C 2.3 S. Valcke 99/04/30 added: printing levels 00059 C 2.4 S. Legutke 00/07/12 conservation calculation 00060 C 3.0 V. Gayler 07/01/23 added GLBPOS conserv. method 00061 C 00062 C %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 00063 C 00064 C* ---------------------------- Include files --------------------------- 00065 C 00066 USE mod_unit 00067 USE mod_printing 00068 C 00069 C* ---------------------------- Argument declarations ------------------- 00070 C 00071 REAL (kind=ip_realwp_p) pflda(ksiza), pfldb(ksizb) 00072 REAL (kind=ip_realwp_p) psurfa(ksiza), psurfb(ksizb) 00073 REAL (kind=ip_realwp_p) psgra(ksiza), psgrb(ksizb) 00074 INTEGER (kind=ip_intwp_p) kmska(ksiza), kmskb(ksizb) 00075 CHARACTER*8 cdmet 00076 C 00077 C* ---------------------------- Poema verses ---------------------------- 00078 C 00079 C %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 00080 C 00081 C* 1. Initialization 00082 C -------------- 00083 C 00084 IF (nlogprt .GE. 2) THEN 00085 WRITE (UNIT = nulou,FMT = *) ' ' 00086 WRITE (UNIT = nulou,FMT = *) ' ' 00087 WRITE (UNIT = nulou,FMT = *) 00088 $ ' ROUTINE conserv - Level 3' 00089 WRITE (UNIT = nulou,FMT = *) 00090 $ ' *************** *******' 00091 WRITE (UNIT = nulou,FMT = *) ' ' 00092 WRITE (UNIT = nulou,FMT = *) 00093 $ ' Global or local flux conservation' 00094 WRITE (UNIT = nulou,FMT = *) ' ' 00095 CALL FLUSH(nulou) 00096 ENDIF 00097 zrad = 6371229. 00098 zrad2 = zrad * zrad 00099 zradi = 1./zrad2 00100 00101 C 00102 C 00103 C* 2. Calculate mesh surfaces 00104 C ----------------------- 00105 C 00106 DO 210 ji = 1, ksiza 00107 psgra(ji) = psurfa(ji) * zradi 00108 210 CONTINUE 00109 DO 220 ji = 1, ksizb 00110 psgrb(ji) = psurfb(ji) * zradi 00111 220 CONTINUE 00112 C 00113 C 00114 C* 3. Lagrange multiplier for global conservation 00115 C ------------------------------------------- 00116 C 00117 IF (cdmet .EQ. 'GLOBAL' .OR. cdmet .EQ. 'GLBPOS' 00118 $ .OR. cdmet .EQ. 'BASBAL' .OR. cdmet .EQ. 'BASPOS' ) THEN 00119 zflxa = 0. 00120 ztsqi = 0. 00121 ztsqs = 0. 00122 C 00123 C* Sum up flux on target grid 00124 C 00125 DO 310 ji = 1, ksiza 00126 IF (kmska(ji) .eq. 0) THEN 00127 zflxa = zflxa + psgra(ji) * pflda(ji) 00128 CSL ztsqi = ztsqi + psgra(ji) * psgra(ji) 00129 ztsqi = ztsqi + psgra(ji) 00130 ENDIF 00131 310 CONTINUE 00132 C 00133 C* Sum up flux on source grid 00134 C 00135 zflxb = 0. 00136 DO 320 ji = 1, ksizb 00137 IF (kmskb(ji) .eq. 0) THEN 00138 zflxb = zflxb + psgrb(ji) * pfldb(ji) 00139 ztsqs = ztsqs + psgrb(ji) 00140 ENDIF 00141 320 CONTINUE 00142 C 00143 IF (cdmet .EQ. 'GLOBAL') THEN 00144 C 00145 C* Get global correction 00146 C 00147 zlagr = (zflxa - zflxb) / ztsqi 00148 C 00149 C* Constrained solution: the error is uniformly shared between sea points 00150 C 00151 DO 330 ji = 1, ksiza 00152 IF (kmska(ji) .EQ. 0) THEN 00153 CSL pflda(ji) = pflda(ji) - zlagr * psgra(ji) 00154 pflda(ji) = pflda(ji) - zlagr 00155 ENDIF 00156 330 CONTINUE 00157 C 00158 ELSE IF (cdmet .EQ. 'GLBPOS') THEN 00159 C 00160 C* Get global correction that does not change field sign 00161 C 00162 IF (zflxa .EQ. 0. .AND. zflxb .NE. 0.) THEN 00163 CALL HALTE('STOP in conserv: zflxra = 0') 00164 ELSE IF (zflxa .NE. 0.) THEN 00165 zlagr = zflxb / zflxa 00166 00167 DO 333 ji = 1, ksiza 00168 IF (kmska(ji) .EQ. 0) THEN 00169 pflda(ji) = pflda(ji) * zlagr 00170 ENDIF 00171 333 CONTINUE 00172 ENDIF 00173 C 00174 ELSE IF (cdmet .EQ. 'BASBAL') THEN 00175 C 00176 C* Get global correction 00177 C 00178 zlagr = (zflxa - zflxb*(ztsqi/ztsqs)) / ztsqi 00179 C 00180 C* Constrained solution: the error is uniformly shared between sea points 00181 C 00182 DO 336 ji = 1, ksiza 00183 IF (kmska(ji) .EQ. 0) THEN 00184 pflda(ji) = pflda(ji) - zlagr 00185 ENDIF 00186 336 CONTINUE 00187 C 00188 ELSE IF (cdmet .EQ. 'BASPOS') THEN 00189 C 00190 C* Get global correction that does not change field sign 00191 C 00192 IF (zflxa .EQ. 0. .AND. zflxb .NE. 0.) THEN 00193 CALL HALTE('STOP in conserv: zflxra = 0') 00194 ELSE IF (zflxa .NE. 0.) THEN 00195 zlagr = (zflxb / zflxa)*(ztsqi/ztsqs) 00196 00197 DO 338 ji = 1, ksiza 00198 IF (kmska(ji) .EQ. 0) THEN 00199 pflda(ji) = pflda(ji) * zlagr 00200 ENDIF 00201 338 CONTINUE 00202 ENDIF 00203 ENDIF 00204 C 00205 C* Printing test 00206 C 00207 zflxn = 0. 00208 DO 340 ji = 1, ksiza 00209 IF (kmska(ji) .eq. 0) THEN 00210 zflxn = zflxn + psgra(ji) * pflda(ji) 00211 ENDIF 00212 340 CONTINUE 00213 IF (nlogprt .GE. 2) THEN 00214 WRITE (UNIT = nulou,FMT = *) 00215 $ ' Printing check for flux conservation ' 00216 WRITE (UNIT = nulou,FMT = *) ' ' 00217 WRITE (UNIT = nulou,FMT = *) 00218 $ ' Total surface on source grid ZTsQS = ',ztsqs 00219 WRITE (UNIT = nulou,FMT = *) 00220 $ ' Total surface on target grid ZTsQI = ',ztsqi 00221 WRITE (UNIT = nulou,FMT = *) 00222 $ ' Total flux on source grid ZFLXB = ',zflxb 00223 WRITE (UNIT = nulou,FMT = *) 00224 $ ' Total flux on target grid ZFLXA = ',zflxa 00225 WRITE (UNIT = nulou,FMT = *) 00226 $ ' Difference ZLAGR = ',zlagr 00227 WRITE (UNIT = nulou,FMT = *) 00228 $ ' Idem after conservation ZFLXN = ',zflxn 00229 WRITE (UNIT = nulou,FMT = *) ' ' 00230 ENDIF 00231 ENDIF 00232 C 00233 C 00234 C* 4. End of routine 00235 C -------------- 00236 C 00237 IF (nlogprt .GE. 2) THEN 00238 WRITE (UNIT = nulou,FMT = *) ' ' 00239 WRITE (UNIT = nulou,FMT = *) 00240 $ ' --------- End of ROUTINE conserv ---------' 00241 CALL FLUSH (nulou) 00242 ENDIF 00243 RETURN 00244 END
1.7.4 
