1313! ! subroutine const_diffusive_profile_fixS(is_scalar, is_source, &
1414! ! & sph_rj, r_2nd, sc_prop, sph_bc, bcs_S, fdm2_center, &
1515! ! & band_s00_poisson, rj_fld, reftemp_rj, reftemp_local)
16- ! ! subroutine gradient_of_radial_reference(sph_rj, r_2nd, &
17- ! ! & sph_bc, bcs_S , fdm2_center, reftemp_r, refgrad_r)
16+ ! ! subroutine gradient_of_radial_reference(sph_rj, sph_bc, bcs_S, &
17+ ! ! & r_2nd , fdm2_center, reftemp_r, refgrad_r)
1818! ! integer(kind = kint), intent(in) :: is_scalar, is_source
1919! ! type(sph_rj_grid), intent(in) :: sph_rj
2020! ! type(fdm_matrices), intent(in) :: r_2nd
@@ -71,11 +71,7 @@ subroutine s_const_diffusive_profile(sph_rj, r_2nd, &
7171!
7272 use calypso_mpi
7373 use calypso_mpi_real
74- use const_sph_radial_grad
7574 use fill_scalar_field
76- use select_exp_scalar_ICB
77- use select_exp_scalar_CMB
78- ! use cal_sol_reftemp_BiCGSTAB
7975!
8076 type (sph_rj_grid), intent (in ) :: sph_rj
8177 type (fdm_matrices), intent (in ) :: r_2nd
@@ -93,144 +89,100 @@ subroutine s_const_diffusive_profile(sph_rj, r_2nd, &
9389 & :: grad_local(0 :sph_rj% nidx_rj(1 ))
9490!
9591 integer (kind = kint_gl) :: num64
96- ! integer :: k
9792!
9893 if (sph_rj% idx_rj_degree_zero .gt. 0 ) then
99- ! $omp parallel workshare
100- ref_local(0 :sph_rj% nidx_rj(1 )) = ref_local(0 :sph_rj% nidx_rj(1 )) &
101- & * (sc_prop% coef_source / sc_prop% coef_diffuse)
102- ! $omp end parallel workshare
103- !
104- call set_ICB_scalar_boundary_1d &
105- & (sph_rj, sph_bc, bcs_S% ICB_Sspec, ref_local(0 ))
106- call set_CMB_scalar_boundary_1d &
107- & (sph_rj, sph_bc, bcs_S% CMB_Sspec, ref_local(0 ))
108- !
109- ! do k = 0, sph_rj%nidx_rj(1)
110- ! write(*,*) k, 'RHS', ref_local(k)
111- ! end do
112- !
113- call lubksb_3band_ctr(band_s00_poisson, ref_local(0 ))
114- ! call s_cal_sol_reftemp_BiCGSTAB &
115- ! & (band_s00_poisson, ref_local(0))
116- !
117- ! do k = 0, sph_rj%nidx_rj(1)
118- ! write(*,*) k, 'Solution', ref_local(k)
119- ! end do
120- !
94+ call cal_diffusive_profile &
95+ & (sph_rj, sc_prop, sph_bc, bcs_S, r_2nd, fdm2_center, &
96+ & band_s00_poisson, ref_local)
12197 call fill_scalar_1d_external(sph_bc, sph_rj% inod_rj_center, &
122- & sph_rj% nidx_rj(1 ), ref_local(0 ))
123- !
124- call cal_sph_nod_gradient_1d(sph_bc% kr_in, sph_bc% kr_out, &
125- & sph_rj% nidx_rj(1 ), r_2nd% fdm(1 )% dmat, &
126- & ref_local(0 ), grad_local(0 ))
98+ & sph_rj% nidx_rj(1 ), ref_local)
12799!
128- call sel_ICB_radial_grad_1d_scalar &
129- & (sph_rj, sph_bc, bcs_S% ICB_Sspec, fdm2_center, &
130- & ref_local, grad_local)
131- call sel_CMB_radial_grad_1d_scalar &
132- & (sph_rj, sph_bc, bcs_S% CMB_Sspec, ref_local, grad_local)
100+ call gradient_of_radial_reference(sph_rj, sph_bc, bcs_S, &
101+ & r_2nd, fdm2_center, ref_local, grad_local)
133102 end if
134103!
135104! $omp parallel workshare
136105 reftemp_r(0 :sph_rj% nidx_rj(1 )) = 0.0d0
137- refgrad_r(0 :sph_rj% nidx_rj(1 )) = 0.0d0
138106! $omp end parallel workshare
139107!
140108 num64 = sph_rj% nidx_rj(1 ) + 1
141109 call calypso_mpi_allreduce_real(ref_local(0 ), reftemp_r, &
142110 & num64, MPI_SUM)
111+ !
112+ !
113+ ! $omp parallel workshare
114+ refgrad_r(0 :sph_rj% nidx_rj(1 )) = 0.0d0
115+ ! $omp end parallel workshare
116+ !
117+ num64 = sph_rj% nidx_rj(1 ) + 1
143118 call calypso_mpi_allreduce_real(grad_local(0 ), refgrad_r, &
144119 & num64, MPI_SUM)
145120!
146121 end subroutine s_const_diffusive_profile
147122!
148123! -----------------------------------------------------------------------
149124!
150- subroutine const_diffusive_profile_fixS ( is_scalar , is_source , &
151- & sph_rj , r_2nd , sc_prop , sph_bc , bcs_S , fdm2_center , &
152- & band_s00_poisson , rj_fld , reftemp_rj , reftemp_local )
125+ subroutine cal_diffusive_profile &
126+ & ( sph_rj, sc_prop, sph_bc, bcs_S, r_2nd , fdm2_center, &
127+ & band_s00_poisson, ref_local )
153128!
154- use calypso_mpi
155- use calypso_mpi_real
156129 use const_sph_radial_grad
157130 use fill_scalar_field
158131 use select_exp_scalar_ICB
159132 use select_exp_scalar_CMB
133+ ! use cal_sol_reftemp_BiCGSTAB
160134!
161- integer (kind = kint), intent (in ) :: is_scalar, is_source
162135 type (sph_rj_grid), intent (in ) :: sph_rj
163136 type (fdm_matrices), intent (in ) :: r_2nd
164137 type (scalar_property), intent (in ) :: sc_prop
165138 type (sph_boundary_type), intent (in ) :: sph_bc
166139 type (sph_scalar_boundary_data), intent (in ) :: bcs_S
167140 type (fdm2_center_mat), intent (in ) :: fdm2_center
168- type (phys_data), intent (in ) :: rj_fld
169141 type (band_matrix_type), intent (in ) :: band_s00_poisson
170142!
171143 real (kind = kreal), intent (inout ) &
172- & :: reftemp_rj(0 :sph_rj% nidx_rj(1 ),0 :1 )
173- real (kind = kreal), intent (inout ) &
174- & :: reftemp_local(0 :sph_rj% nidx_rj(1 ),0 :1 )
144+ & :: ref_local(0 :sph_rj% nidx_rj(1 ))
175145!
176- integer (kind = kint) :: inod
177- integer (kind = kint_gl) :: num64
146+ ! real(kind = kreal), allocatable :: ref_CG(:)
178147!
148+ ! integer :: k
179149!
180- if (sph_rj% idx_rj_degree_zero .gt. 0 ) then
181- if (is_source .gt. 0 ) then
182- call copy_degree0_comps_to_sol(sph_rj% nidx_rj(2 ), &
183- & sph_rj% inod_rj_center, sph_rj% idx_rj_degree_zero, &
184- & rj_fld% n_point, rj_fld% d_fld(1 ,is_source), &
185- & sph_rj% nidx_rj(1 ), reftemp_local(0 ,0 ))
186- ! $omp parallel workshare
187- reftemp_local(0 :sph_rj% nidx_rj(1 ),0 ) &
188- & = (sc_prop% coef_source / sc_prop% coef_diffuse) &
189- & * reftemp_local(0 :sph_rj% nidx_rj(1 ),0 )
190- ! $omp end parallel workshare
191- else
192150! $omp parallel workshare
193- reftemp_local(0 :sph_rj% nidx_rj(1 ),0 ) = 0.0d0
151+ ref_local(0 :sph_rj% nidx_rj(1 )) = ref_local(0 :sph_rj% nidx_rj(1 )) &
152+ & * (sc_prop% coef_source / sc_prop% coef_diffuse)
194153! $omp end parallel workshare
195- end if
196154!
197- call set_ICB_scalar_boundary_1d &
198- & (sph_rj, sph_bc, bcs_S% ICB_Sspec, reftemp_local(0 ,0 ))
199- !
200- inod = sph_rj% idx_rj_degree_zero &
201- & + (sph_bc% kr_out-1 ) * sph_rj% nidx_rj(2 )
202- reftemp_local(sph_bc% kr_out,0 ) = rj_fld% d_fld(inod,is_scalar)
155+ call set_ICB_scalar_boundary_1d(sph_rj, sph_bc, bcs_S% ICB_Sspec, &
156+ & ref_local(0 ))
157+ call set_CMB_scalar_boundary_1d(sph_rj, sph_bc, bcs_S% CMB_Sspec, &
158+ & ref_local(0 ))
203159!
204- call lubksb_3band_ctr(band_s00_poisson, reftemp_local(0 ,0 ))
205- call fill_scalar_1d_external(sph_bc, sph_rj% inod_rj_center, &
206- & sph_rj% nidx_rj(1 ), reftemp_local(0 ,0 ))
160+ call lubksb_3band_ctr(band_s00_poisson, ref_local(0 ))
207161!
208- call cal_sph_nod_gradient_1d(sph_bc% kr_in, sph_bc% kr_out, &
209- & sph_rj% nidx_rj(1 ), r_2nd% fdm(1 )% dmat, &
210- & reftemp_local(0 ,0 ), reftemp_local(0 ,1 ))
162+ ! allocate(ref_CG(0:sph_rj%nidx_rj(1)))
163+ ! !$omp parallel workshare
164+ ! ref_CG(0:sph_rj%nidx_rj(1)) = ref_local(0:sph_rj%nidx_rj(1))
165+ ! !$omp end parallel workshare
166+ ! call s_cal_sol_reftemp_BiCGSTAB(band_s00_poisson, &
167+ ! & ref_CG(0))
211168!
212- call fix_ICB_radial_grad_1d_scalar(sph_rj, sph_bc, fdm2_center, &
213- & reftemp_local(0 ,0 ), reftemp_local(0 ,1 ))
214- call fix_CMB_radial_grad_1d_scalar(sph_rj, sph_bc, &
215- & reftemp_local(0 ,0 ), reftemp_local(0 ,1 ))
216- end if
169+ ! do k = 0, sph_rj%nidx_rj(1)
170+ ! write(*,*) k, 'RHS', ref_local(k)
171+ ! end do
217172!
218- ! $omp parallel workshare
219- reftemp_rj(0 :sph_rj% nidx_rj(1 ),0 :1 ) = 0.0d0
220- ! $omp end parallel workshare
221- num64 = 2 * (sph_rj% nidx_rj(1 ) + 1 )
222- call calypso_mpi_allreduce_real(reftemp_local, reftemp_rj, &
223- & num64, MPI_SUM)
173+ ! do k = 0, sph_rj%nidx_rj(1)
174+ ! write(*,*) k, 'Solution', ref_local(k), ref_CG(k)
175+ ! end do
176+ ! deallocate(ref_CG)
224177!
225- end subroutine const_diffusive_profile_fixS
178+ end subroutine cal_diffusive_profile
226179!
227180! -----------------------------------------------------------------------
228181!
229- subroutine gradient_of_radial_reference (sph_rj , r_2nd , &
230- & sph_bc , bcs_S , fdm2_center , reftemp_r , refgrad_r )
182+ subroutine gradient_of_radial_reference (sph_rj , sph_bc , bcs_S , &
183+ & r_2nd , fdm2_center , reftemp_r , refgrad_r )
231184!
232185 use const_sph_radial_grad
233- use fill_scalar_field
234186 use select_exp_scalar_ICB
235187 use select_exp_scalar_CMB
236188!
@@ -239,13 +191,10 @@ subroutine gradient_of_radial_reference(sph_rj, r_2nd, &
239191 type (sph_boundary_type), intent (in ) :: sph_bc
240192 type (sph_scalar_boundary_data), intent (in ) :: bcs_S
241193 type (fdm2_center_mat), intent (in ) :: fdm2_center
194+ real (kind= kreal), intent (in ) :: reftemp_r(0 :sph_rj% nidx_rj(1 ))
242195!
243- real (kind= kreal), intent (inout ) :: reftemp_r(0 :sph_rj% nidx_rj(1 ))
244196 real (kind= kreal), intent (inout ) :: refgrad_r(0 :sph_rj% nidx_rj(1 ))
245197!
246- !
247- call fill_scalar_1d_external(sph_bc, sph_rj% inod_rj_center, &
248- & sph_rj% nidx_rj(1 ), reftemp_r(0 ))
249198!
250199 call cal_sph_nod_gradient_1d(sph_bc% kr_in, sph_bc% kr_out, &
251200 & sph_rj% nidx_rj(1 ), r_2nd% fdm(1 )% dmat, &
@@ -277,11 +226,11 @@ subroutine cal_reference_source(sph_rj, sc_prop, &
277226!
278227!
279228 ref_local(0 :sph_rj% nidx_rj(1 )) = 0.0d0
280- k = 1
229+ k = 0
281230 ref_local(k) = band_s00_poisson% mat(2 ,k ) * ref_source(k) &
282231 & + band_s00_poisson% mat(1 ,k+1 ) * ref_source(k+1 )
283232! $omp parallel do
284- do k = 2 , sph_rj% nidx_rj(1 ) - 1
233+ do k = 1 , sph_rj% nidx_rj(1 ) - 1
285234 ref_local(k) = band_s00_poisson% mat(3 ,k-1 ) * ref_source(k-1 ) &
286235 & + band_s00_poisson% mat(2 ,k ) * ref_source(k) &
287236 & + band_s00_poisson% mat(1 ,k+1 ) * ref_source(k+1 )
@@ -299,5 +248,84 @@ subroutine cal_reference_source(sph_rj, sc_prop, &
299248 end subroutine cal_reference_source
300249!
301250! -----------------------------------------------------------------------
251+ !
252+ subroutine const_diffusive_profile_fixS (is_scalar , is_source , &
253+ & sph_rj , r_2nd , sc_prop , sph_bc , bcs_S , fdm2_center , &
254+ & band_s00_poisson , rj_fld , reftemp_rj , reftemp_local )
255+ !
256+ use calypso_mpi
257+ use calypso_mpi_real
258+ use const_sph_radial_grad
259+ use fill_scalar_field
260+ use select_exp_scalar_ICB
261+ use select_exp_scalar_CMB
262+ !
263+ integer (kind = kint), intent (in ) :: is_scalar, is_source
264+ type (sph_rj_grid), intent (in ) :: sph_rj
265+ type (fdm_matrices), intent (in ) :: r_2nd
266+ type (scalar_property), intent (in ) :: sc_prop
267+ type (sph_boundary_type), intent (in ) :: sph_bc
268+ type (sph_scalar_boundary_data), intent (in ) :: bcs_S
269+ type (fdm2_center_mat), intent (in ) :: fdm2_center
270+ type (phys_data), intent (in ) :: rj_fld
271+ type (band_matrix_type), intent (in ) :: band_s00_poisson
272+ !
273+ real (kind = kreal), intent (inout ) &
274+ & :: reftemp_rj(0 :sph_rj% nidx_rj(1 ),0 :1 )
275+ real (kind = kreal), intent (inout ) &
276+ & :: reftemp_local(0 :sph_rj% nidx_rj(1 ),0 :1 )
277+ !
278+ integer (kind = kint) :: inod
279+ integer (kind = kint_gl) :: num64
280+ !
281+ !
282+ if (sph_rj% idx_rj_degree_zero .gt. 0 ) then
283+ if (is_source .gt. 0 ) then
284+ call copy_degree0_comps_to_sol(sph_rj% nidx_rj(2 ), &
285+ & sph_rj% inod_rj_center, sph_rj% idx_rj_degree_zero, &
286+ & rj_fld% n_point, rj_fld% d_fld(1 ,is_source), &
287+ & sph_rj% nidx_rj(1 ), reftemp_local(0 ,0 ))
288+ ! $omp parallel workshare
289+ reftemp_local(0 :sph_rj% nidx_rj(1 ),0 ) &
290+ & = (sc_prop% coef_source / sc_prop% coef_diffuse) &
291+ & * reftemp_local(0 :sph_rj% nidx_rj(1 ),0 )
292+ ! $omp end parallel workshare
293+ else
294+ ! $omp parallel workshare
295+ reftemp_local(0 :sph_rj% nidx_rj(1 ),0 ) = 0.0d0
296+ ! $omp end parallel workshare
297+ end if
298+ !
299+ call set_ICB_scalar_boundary_1d &
300+ & (sph_rj, sph_bc, bcs_S% ICB_Sspec, reftemp_local(0 ,0 ))
301+ !
302+ inod = sph_rj% idx_rj_degree_zero &
303+ & + (sph_bc% kr_out-1 ) * sph_rj% nidx_rj(2 )
304+ reftemp_local(sph_bc% kr_out,0 ) = rj_fld% d_fld(inod,is_scalar)
305+ !
306+ call lubksb_3band_ctr(band_s00_poisson, reftemp_local(0 ,0 ))
307+ call fill_scalar_1d_external(sph_bc, sph_rj% inod_rj_center, &
308+ & sph_rj% nidx_rj(1 ), reftemp_local(0 ,0 ))
309+ !
310+ call cal_sph_nod_gradient_1d(sph_bc% kr_in, sph_bc% kr_out, &
311+ & sph_rj% nidx_rj(1 ), r_2nd% fdm(1 )% dmat, &
312+ & reftemp_local(0 ,0 ), reftemp_local(0 ,1 ))
313+ !
314+ call fix_ICB_radial_grad_1d_scalar(sph_rj, sph_bc, fdm2_center, &
315+ & reftemp_local(0 ,0 ), reftemp_local(0 ,1 ))
316+ call fix_CMB_radial_grad_1d_scalar(sph_rj, sph_bc, &
317+ & reftemp_local(0 ,0 ), reftemp_local(0 ,1 ))
318+ end if
319+ !
320+ ! $omp parallel workshare
321+ reftemp_rj(0 :sph_rj% nidx_rj(1 ),0 :1 ) = 0.0d0
322+ ! $omp end parallel workshare
323+ num64 = 2 * (sph_rj% nidx_rj(1 ) + 1 )
324+ call calypso_mpi_allreduce_real(reftemp_local, reftemp_rj, &
325+ & num64, MPI_SUM)
326+ !
327+ end subroutine const_diffusive_profile_fixS
328+ !
329+ ! -----------------------------------------------------------------------
302330!
303331 end module const_diffusive_profile
0 commit comments