Newer
Older
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
subroutine read_restart(params,istep,iter,current_time,osolve,ov,vo,surface,cl,&
zi,nest)
use definitions
!use mpi
implicit none
include 'mpif.h'
! Passed variables
type (parameters) :: params
integer :: istep,iter
double precision :: current_time
type (octreesolve) :: osolve
type (octreev) :: ov
type (void) :: vo
type (sheet) :: surface(params%ns)
type (cloud) :: cl
type (ziso) :: zi
type (nest_info) :: nest
! Internal variables
double precision :: current_time
integer :: err;i,is,k
osolve%noctree=params%noctreemax
allocate (osolve%octree(osolve%noctree),stat=err) ; if (err.ne.0) call stop_run ('Error alloc osolve%octree in read_restart$')
open (9,file=trim(params%restartfile),status='old',form='unformatted')
read (9) osolve%octree(3), &
osolve%nnode, &
osolve%nleaves, &
osolve%nface, &
osolve%nlsf, &
cl%np, &
current_time
! Massive allocation
allocate (osolve%x(osolve%nnode),stat=err) ; if (err.ne.0) call stop_run ('Error alloc osolve%x in read_restart$')
allocate (osolve%y(osolve%nnode),stat=err) ; if (err.ne.0) call stop_run ('Error alloc osolve%y in read_restart$')
allocate (osolve%z(osolve%nnode),stat=err) ; if (err.ne.0) call stop_run ('Error alloc osolve%z in read_restart$')
allocate (osolve%unode(osolve%nnode),stat=err) ; if (err.ne.0) call stop_run ('Error alloc osolve%unode in read_restart$')
allocate (osolve%vnode(osolve%nnode),stat=err) ; if (err.ne.0) call stop_run ('Error alloc osolve%vnode in read_restart$')
allocate (osolve%wnode(osolve%nnode),stat=err) ; if (err.ne.0) call stop_run ('Error alloc osolve%wnode in read_restart$')
allocate (osolve%wnodeiso(osolve%nnode),stat=err) ; if (err.ne.0) call stop_run ('Error alloc osolve%wnodeiso in read_restart$')
allocate (osolve%lsf(osolve%nnode,osolve%nlsf),stat=err) ; if (err.ne.0) call stop_run ('Error alloc osolve%lsf in read_restart$')
allocate (osolve%temp(osolve%nnode),stat=err) ; if (err.ne.0) call stop_run ('Error alloc osolve%temp in read_restart$')
allocate (osolve%strain(osolve%nnode),stat=err) ; if (err.ne.0) call stop_run ('Error alloc osolve%strain in read_restart$')
allocate (osolve%kfix(osolve%nnode*3),stat=err) ; if (err.ne.0) call stop_run ('Error alloc osolve%kfix in read_restart$')
allocate (osolve%kfixt(osolve%nnode),stat=err) ; if (err.ne.0) call stop_run ('Error alloc osolve%kfixt in read_restart$')
allocate (osolve%icon(8,osolve%nleaves),stat=err) ; if (err.ne.0) call stop_run ('Error alloc osolve%icon in read_restart$')
allocate (osolve%pressure(osolve%nleaves),stat=err) ; if (err.ne.0) call stop_run ('Error alloc osolve%pressure in read_restart$')
allocate (osolve%spressure(osolve%nleaves),stat=err) ; if (err.ne.0) call stop_run ('Error alloc osolve%spressure in read_restart$')
allocate (osolve%crit(osolve%nleaves),stat=err) ; if (err.ne.0) call stop_run ('Error alloc osolve%crit in read_restart$')
allocate (osolve%e2d(osolve%nleaves),stat=err) ; if (err.ne.0) call stop_run ('Error alloc osolve%e2d in read_restart$')
allocate (osolve%eviscosity(osolve%nleaves),stat=err) ; if (err.ne.0) call stop_run ('Error alloc osolve%eviscosity in read_restart$')
allocate (osolve%is_plastic(osolve%nleaves),stat=err) ; if (err.ne.0) call stop_run ('Error alloc osolve%is_plastic in read_restart$')
allocate (osolve%dilatr(osolve%nleaves),stat=err) ; if (err.ne.0) call stop_run ('Error alloc osolve%dilatr in read_restart$')
allocate (osolve%matnum(osolve%nleaves),stat=err) ; if (err.ne.0) call stop_run ('Error alloc osolve%matnum in read_restart$')
allocate (osolve%iface(9,osolve%nface),stat=err) ; if (err.ne.0) call stop_run ('Error alloc osolve%iface in read_restart$')
allocate (vo%node(osolve%nnode),stat=threadinfo%err) ; if (err.ne.0) call stop_run ('Error alloc vo%node in read_restart$')
allocate (vo%leaf(osolve%nnode),stat=threadinfo%err) ; if (err.ne.0) call stop_run ('Error alloc vo%leaf in read_restart$')
allocate (vo%ftr(osolve%nnode),stat=threadinfo%err) ; if (err.ne.0) call stop_run ('Error alloc vo%ftr in read_restart$')
allocate (vo%rtf(osolve%nnode),stat=threadinfo%err) ; if (err.ne.0) call stop_run ('Error alloc vo%rtf in read_restart$')
allocate (vo%influid(osolve%nnode),stat=threadinfo%err) ; if (err.ne.0) call stop_run ('Error alloc vo%influid in read_restart$')
allocate (vo%face(osolve%nface),stat=threadinfo%err) ; if (err.ne.0) call stop_run ('Error alloc vo%face in read_restart$')
allocate (ov%temporary_nodal_pressure(ov%nnode),stat=err) ; if (err.ne.0) call stop_run ('Error alloc ov%temp_nodal_pressure in define_ov$')
allocate (ov%whole_leaf_in_fluid(ov%nleaves),stat=err) ; if (err.ne.0) call stop_run ('Error alloc ov%whole_leaf_in_fluid in define_ov$')
! Read info for octree solve nodes (x,y,z,u,v,w,lsf,temp)
read (9) (osolve%x(i), &
osolve%y(i), &
osolve%z(i), &
osolve%u(i), &
osolve%v(i), &
osolve%w(i), &
osolve%wiso(i), &
(osolve%lsf(i,j),j=1,osolve%nlsf), &
osolve%temp(i), &
ov%temporary_nodal_pressure(i), &
osolve%strain(i), &
osolve%kfix((i-1)*3+1), &
osolve%kfix((i-1)*3+2), &
osolve%kfix((i-1)*3+3), &
osolve%kfixt(i), &
i=1,osolve%nnode)
! Read icon array
read (9) ((osolve%icon(k,i),k=1,8), &
osolve%pressure(i), &
osolve%spressure(i), &
osolve%crit(i), &
osolve%e2d(i), &
osolve%eviscosity(i), &
osolve%is_plastic(i), &
osolve%dilatr(i), &
osolve%matnum(i), &
ov%whole_leaf_in_fluid(i), &
i=1,osolve%nleaves)
! Read octree information
read (9) (osolve%octree(i),i=1,osolve%octree(3))
! Read bad face information
read (9) ((osolve%iface(k,i),k=1,9),i=1,osolve%nface)
! Read void information
read (9) (vo%node(i), &
vo%leaf(i), &
vo%ftr(i), &
vo%rtf(i), &
vo%influid(i), &
i=1,osolve%nnode)
! Read bad faces information
read (9) (vo%face(i),i=1,osolve%nface)
! Read surface information (r,s,x,y,z,xn,yn,zn)
do is=1,osolve%nlsf
read (9) surface(is)%nsurface, &
surface(is)%activation_time, &
surface(is)%nt
allocate (surface(is)%r(surface(is)%nsurface),stat=err) ; if (err.ne.0) call stop_run ('Error alloc surface(is)%r in read_restart$')
allocate (surface(is)%s(surface(is)%nsurface),stat=err) ; if (err.ne.0) call stop_run ('Error alloc surface(is)%s in read_restart$')
allocate (surface(is)%x(surface(is)%nsurface),stat=err) ; if (err.ne.0) call stop_run ('Error alloc surface(is)%x in read_restart$')
allocate (surface(is)%y(surface(is)%nsurface),stat=err) ; if (err.ne.0) call stop_run ('Error alloc surface(is)%y in read_restart$')
allocate (surface(is)%z(surface(is)%nsurface),stat=err) ; if (err.ne.0) call stop_run ('Error alloc surface(is)%z in read_restart$')
allocate (surface(is)%xn(surface(is)%nsurface),stat=err) ; if (err.ne.0) call stop_run ('Error alloc surface(is)%xn in read_restart$')
allocate (surface(is)%yn(surface(is)%nsurface),stat=err) ; if (err.ne.0) call stop_run ('Error alloc surface(is)%yn in read_restart$')
allocate (surface(is)%zn(surface(is)%nsurface),stat=err) ; if (err.ne.0) call stop_run ('Error alloc surface(is)%zn in read_restart$')
allocate (surface(is)%u(surface(is)%nsurface),stat=err) ; if (err.ne.0) call stop_run ('Error alloc surface(is)%u in read_restart$')
allocate (surface(is)%v(surface(is)%nsurface),stat=err) ; if (err.ne.0) call stop_run ('Error alloc surface(is)%v in read_restart$')
allocate (surface(is)%w(surface(is)%nsurface),stat=err) ; if (err.ne.0) call stop_run ('Error alloc surface(is)%w in read_restart$')
allocate (surface(is)%icon(3,surface(is)%nt),stat=err) ; if (err.ne.0) call stop_run ('Error alloc surface(is)%icon in read_restart$')
read (9) (surface(is)%r(i), &
surface(is)%s(i), &
surface(is)%x(i), &
surface(is)%y(i), &
surface(is)%z(i), &
surface(is)%xn(i), &
surface(is)%yn(i), &
surface(is)%zn(i), &
surface(is)%u(i), &
surface(is)%v(i), &
surface(is)%w(i), &
i=1,surface(is)%nsurface)
read (9) ((surface(is)%icon(k,i),k=1,3),i=1,surface(is)%nt)
! correct surface for vertical exageration
surface(is)%z=surface(is)%z/params%vex
surface(is)%zn=surface(is)%zn/params%vex
enddo
! Read cloud information
read (9) (cl%x(i), &
cl%y(i), &
cl%z(i), &
cl%x0(i), &
cl%y0(i), &
cl%z0(i), &
cl%strain(i), &
cl%lsf0(i), &
cl%temp(i), &
cl%press(i), &
cl%tag(i), &
i=1,cl%np)
! correct for vertical exaggeration
cl%z=cl%z/params%vex
cl%z0=cl%z0/params%vex
! Read isostasy basal displacement array
read (9) 2**params%levelmax_oct
allocate (zi%zisodisp(2**params%levelmax_oct+1,2**params%levelmax_oct+1), stat=err) ; if (err.ne.0) call stop_run ('Error alloc zi%zisodisp in read_restart$')
if (params%isobc) then
read (9) ((zi%zisodisp(i,j)+surface(osolve%nlsf)%sp01, &
j=1,2**params%levelmax_oct+1),i=1,2**params%levelmax_oct+1)
else
read (9) ((0.d0+surface(osolve%nlsf)%sp01,j=1,2**params%levelmax_oct+1), &
i=1,2**params%levelmax_oct+1)
endif
! Read nested model info
if (params%nest) then
read (9) nest%sselemx,nest%sselemy,nest%sselemz,nest%xminls,nest%yminls,&
nest%zminls
else
nest%sselemx=0.d0
nest%sselemy=0.d0
nest%sselemz=0.d0
nest%xminls=0.d0
nest%yminls=0.d0
nest%zminls=0.d0
endif
close (9)