39 subroutine borehole(field, x_pos, y_pos, dxi, deta, ch_grid, field_val)
47 real(dp),
dimension(0:JMAX,0:IMAX),
intent(in) :: field
48 real(dp),
intent(in) :: x_pos, y_pos, dxi, deta
49 character (len=*),
intent(in) :: ch_grid
51 real(dp),
intent(out) :: field_val
53 integer(i4b) :: i_kl, i_gr, j_kl, j_gr
54 real(dp) :: real_i, real_j, dist1, dist2, dist3, dist4, &
55 weight1, weight2, weight3, weight4
59 real_i = (x_pos-xi(0)) /dxi
60 real_j = (y_pos-eta(0))/deta
64 i_gr = ceiling(real_i)
65 j_gr = ceiling(real_j)
67 if (i_kl <= 1) stop
' borehole: i_kl <= 1 not allowed!'
68 if (j_kl <= 1) stop
' borehole: j_kl <= 1 not allowed!'
69 if (i_gr >= imax-1) stop
' borehole: i_gr >= IMAX-1 not allowed!'
70 if (j_gr >= jmax-1) stop
' borehole: j_gr >= JMAX-1 not allowed!'
72 dist1 = sqrt( (x_pos-xi(i_kl))**2 + (y_pos-eta(j_kl))**2 )
73 dist2 = sqrt( (x_pos-xi(i_kl))**2 + (y_pos-eta(j_gr))**2 )
74 dist3 = sqrt( (x_pos-xi(i_gr))**2 + (y_pos-eta(j_kl))**2 )
75 dist4 = sqrt( (x_pos-xi(i_gr))**2 + (y_pos-eta(j_gr))**2 )
77 weight1 = 1.0_dp/(dist1+eps)
78 weight2 = 1.0_dp/(dist2+eps)
79 weight3 = 1.0_dp/(dist3+eps)
80 weight4 = 1.0_dp/(dist4+eps)
84 if (ch_grid==
'grid')
then
86 field_val = ( weight1*field(j_kl,i_kl) &
87 +weight2*field(j_gr,i_kl) &
88 +weight3*field(j_kl,i_gr) &
89 +weight4*field(j_gr,i_gr) ) &
90 /(weight1 + weight2 + weight3 + weight4)
92 else if (ch_grid==
'sg_x')
then
95 field_val = ( weight1*0.5_dp*(field(j_kl,i_kl)+field(j_kl,i_kl-1)) &
96 +weight2*0.5_dp*(field(j_gr,i_kl)+field(j_gr,i_kl-1)) &
97 +weight3*0.5_dp*(field(j_kl,i_gr)+field(j_kl,i_gr-1)) &
98 +weight4*0.5_dp*(field(j_gr,i_gr)+field(j_gr,i_gr-1)) ) &
99 /(weight1 + weight2 + weight3 + weight4)
101 else if (ch_grid==
'sg_y')
then
104 field_val = ( weight1*0.5_dp*(field(j_kl,i_kl)+field(j_kl-1,i_kl)) &
105 +weight2*0.5_dp*(field(j_gr,i_kl)+field(j_gr-1,i_kl)) &
106 +weight3*0.5_dp*(field(j_kl,i_gr)+field(j_kl-1,i_gr)) &
107 +weight4*0.5_dp*(field(j_gr,i_gr)+field(j_gr-1,i_gr)) ) &
108 /(weight1 + weight2 + weight3 + weight4)
112 stop
' borehole: Parameter ch_grid has undefined value!'
subroutine borehole(field, x_pos, y_pos, dxi, deta, ch_grid, field_val)
Computation of an arbitrary field quantity for a given borehole position x_pos, y_pos by weighed aver...
Declarations of kind types for SICOPOLIS.
Declarations of global variables for SICOPOLIS (for the ANT domain).
Declarations of global variables for SICOPOLIS.