bilinearintegral_geo.f.html | ![]() |
Source file: bilinearintegral_geo.f | |
Directory: /Users/rjl/git/rjleveque/clawpack-4.6.3/geoclaw/2d/lib | |
Converted: Mon Jan 21 2013 at 20:15:56 using clawcode2html | |
This documentation file will not reflect any later changes in the source file. |
c====================================================================== function bilinearintegral(xim,xip,yjm,yjp,x1,x2,y1,y2,dxx,dyy, & z11,z12,z21,z22) c====================================================================== implicit none * !i/o double precision xim,xip,yjm,yjp,x1,x2,y1,y2,dxx,dyy, & z11,z12,z21,z22 double precision bilinearintegral * !local double precision xlow,ylow,xhi,yhi,area,sumxi,sumeta,a,b,c,d c####################################################################### c c bilinearintegral integrates the bilinear with values z## c over the rectangular region xim <= x <= xip, and c yjm <= y <= yjp c written by David L George c Vancouver, WA April 2010 c####################################################################### * !integrate the portion of the bilinear intersected with the c !rectangular cell analytically. c !find limits of integral (this should already be true?) xlow = max(xim,x1) ylow = max(yjm,y1) xhi = min(xip,x2) yhi = min(yjp,y2) * !find the area of integration area = (yhi-ylow)*(xhi-xlow) sumxi = (xhi + xlow - 2.d0*x1)/dxx sumeta = (yhi + ylow - 2.d0*y1)/dyy * !find coefficients of bilinear a*xi + b*eta + c*xi*eta + d a = z21-z11 b = z12-z11 c = z22-z21-z12+z11 d = z11 bilinearintegral = (0.5d0*(a*sumxi + b*sumeta) & + 0.25d0*c*sumxi*sumeta + d)*area return end function c====================================================================== function bilinearintegral_s(xim,xip,yjm,yjp,x1,x2,y1,y2,dxx,dyy, & z11,z12,z21,z22,Rearth,pi) c====================================================================== implicit none * !i/o double precision xim,xip,yjm,yjp,x1,x2,y1,y2,dxx,dyy, & z11,z12,z21,z22,Rearth,pi double precision bilinearintegral_s * !local double precision xlo,ylo,xhi,yhi,delx,dely,a,b,c,d double precision xdiffhi,xdifflo,ydiffhi,ydifflo,xdiff2,d2r,r2d double precision adsinint,cbsinint c####################################################################### c c bilinearintegral integrates the bilinear with values z## c over the rectangular region xim <= x <= xip, and c yjm <= y <= yjp c integration is actually done on the surface of a sphere c written by David L George c Vancouver, WA April 2010 c####################################################################### * !integrate the portion of the bilinear intersected with the c !rectangular cell analytically. d2r = pi/180.d0 r2d = 180.d0/pi c !find limits of integral (this should already be true?) xlo = max(xim,x1) ylo = max(yjm,y1) xhi = min(xip,x2) yhi = min(yjp,y2) delx = xhi - xlo dely = yhi - ylo * !find terms for the integration xdiffhi = xhi - x1 xdifflo = xlo - x1 ydiffhi = yhi - y1 ydifflo = ylo - y1 xdiff2 = 0.5d0*(xdiffhi**2 - xdifflo**2) cbsinint = (r2d*cos(d2r*yhi) + ydiffhi*sin(d2r*yhi)) & -(r2d*cos(d2r*ylo) + ydifflo*sin(d2r*ylo)) adsinint = r2d*(sin(d2r*yhi) - sin(d2r*ylo)) * !find coefficients of bilinear a*xi + b*eta + c*xi*eta + d a = (z21-z11)/dxx b = (z12-z11)/dyy c = (z22-z21-z12+z11)/(dxx*dyy) d = z11 bilinearintegral_s = ((a*xdiff2 + d*delx)*adsinint & + r2d*(c*xdiff2 + b*delx)*cbsinint)*(Rearth*d2r)**2 return end function