1 |
C $Header: /u/gcmpack/MITgcm/pkg/streamice/streamice_invert_surf_forthick.F,v 1.3 2013/06/21 20:49:50 jmc Exp $ |
2 |
C $Name: $ |
3 |
|
4 |
#include "CPP_OPTIONS.h" |
5 |
#include "STREAMICE_OPTIONS.h" |
6 |
|
7 |
SUBROUTINE STREAMICE_INVERT_SURF_FORTHICK ( |
8 |
O H, |
9 |
I S, |
10 |
I R, |
11 |
I delta, |
12 |
I myThid) |
13 |
|
14 |
! This S/R finds a thickness (H) that gives surf elev S with bed R |
15 |
|
16 |
|
17 |
#include "SIZE.h" |
18 |
#include "GRID.h" |
19 |
#include "SET_GRID.h" |
20 |
#include "EEPARAMS.h" |
21 |
#include "PARAMS.h" |
22 |
#include "STREAMICE.h" |
23 |
|
24 |
#ifdef ALLOW_OPENAD |
25 |
use OAD_tape |
26 |
use OAD_rev |
27 |
use OAD_cp |
28 |
#endif |
29 |
|
30 |
_RL H(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy) |
31 |
_RL S(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy) |
32 |
_RL R(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy) |
33 |
_RL DELTA |
34 |
INTEGER myThid |
35 |
#ifdef ALLOW_OPENAD |
36 |
type(active) :: ETA_GL_STREAMICE |
37 |
type(active) :: X,Y0 |
38 |
type(modeType) :: our_orig_mode |
39 |
#endif |
40 |
|
41 |
|
42 |
#ifdef ALLOW_STREAMICE |
43 |
|
44 |
_RL resid, f, fp, hf, htmp |
45 |
INTEGER i, j, bi, bj, ITER |
46 |
_RL ETA_GL_PRIME_STREAMICE |
47 |
#ifndef ALLOW_OPENAD |
48 |
_RL ETA_GL_STREAMICE |
49 |
EXTERNAL ETA_GL_STREAMICE |
50 |
#endif |
51 |
! EXTERNAL PHI_GL_STREAMICE_PRIME |
52 |
|
53 |
DO bj=myByLo(myThid), myByHi(myThid) |
54 |
DO bi=myBxLo(myThid), myBxHi(myThid) |
55 |
DO j = 1,sNy |
56 |
DO i = 1,sNx |
57 |
IF (STREAMICE_hmask(i,j,bi,bj).eq.1.0) THEN |
58 |
hf = (-1. _d 0) * R(i,j,bi,bj) / |
59 |
& (1. _d 0 - delta) |
60 |
|
61 |
IF (S(i,j,bi,bj) .gt. delta*HF) THEN |
62 |
htmp = S(i,j,bi,bj)-R(i,j,bi,bj) |
63 |
ELSE |
64 |
htmp = S(i,j,bi,bj)/delta |
65 |
ENDIF |
66 |
|
67 |
IF (streamice_smooth_gl_width.gt.0.) THEN |
68 |
|
69 |
RESID=1. _d 0 |
70 |
|
71 |
DO ITER=1,20 |
72 |
IF ((RESID .gt. .005) .and. |
73 |
& ( STREAMICE_hmask(i,j,bi,bj).eq.1.0)) THEN |
74 |
|
75 |
hf = (-1. _d 0) * R(i,j,bi,bj) / |
76 |
& (1. _d 0 - delta) |
77 |
|
78 |
IF (S(i,j,bi,bj) .gt. delta*HF) THEN |
79 |
htmp = S(i,j,bi,bj)-R(i,j,bi,bj) |
80 |
ELSE |
81 |
htmp = S(i,j,bi,bj)/delta |
82 |
ENDIF |
83 |
|
84 |
#ifdef ALLOW_OPENAD |
85 |
|
86 |
our_orig_mode = our_rev_mode |
87 |
our_rev_mode%arg_store=.FALSE. |
88 |
our_rev_mode%arg_restore=.FALSE. |
89 |
our_rev_mode%plain=.TRUE. |
90 |
our_rev_mode%tape=.FALSE. |
91 |
our_rev_mode%adjoint=.FALSE. |
92 |
|
93 |
X%v = htmp-HF |
94 |
Y0%v = delta*HF |
95 |
|
96 |
CALL OpenAD_oad_s_eta_gl_streamice( |
97 |
& X, |
98 |
& delta, |
99 |
& 1. _d 0, |
100 |
& Y0, |
101 |
& streamice_smooth_gl_width, |
102 |
& ETA_GL_STREAMICE) |
103 |
|
104 |
RESID = ETA_GL_STREAMICE%v |
105 |
our_rev_mode = our_orig_mode |
106 |
|
107 |
#else |
108 |
RESID = ETA_GL_STREAMICE ( |
109 |
& htmp-HF, |
110 |
& delta, |
111 |
& 1. _d 0, |
112 |
& delta*HF, |
113 |
& streamice_smooth_gl_width) |
114 |
#endif |
115 |
RESID = RESID - S(i,j,bi,bj) |
116 |
FP = ETA_GL_PRIME_STREAMICE ( |
117 |
& htmp-HF, |
118 |
& delta, |
119 |
& 1. _d 0, |
120 |
& delta*HF, |
121 |
& streamice_smooth_gl_width) |
122 |
Htmp = HTMP-RESID/FP |
123 |
ENDIF |
124 |
ENDDO |
125 |
ENDIF |
126 |
H(i,j,bi,bj) = Htmp |
127 |
ENDIF |
128 |
ENDDO |
129 |
ENDDO |
130 |
ENDDO |
131 |
ENDDO |
132 |
|
133 |
|
134 |
#endif |
135 |
RETURN |
136 |
END |