function[sour,dipl, ier] = rank2d(x,y,dl); % This gives velocity potentials at(x,y), called sour and dipl. % sour is from a unit strength source on dl % dipl is from a unit dipole on dl tol = 0.5e-6; tol1 = 1.0e20; sour = 0.0; dipl = 0.0; x2 = x*x; y2 = y*y; r2 = x2 + y2; d2 = dl * dl; rc = r2/d2; dl2 = 0.5*dl; % Test Input Vlues x, y, dl for overflow if ( (x2 > tol1) | (y2 > tol1) | (d2 > tol1)); ier=1; elseif (dl < tol); ier = 0; elseif (rc > 25); % Far Field Expansions ier = 0; d12 = 1.0/12.0; r2i = 1.0/r2; dr2 = d2 * r2i; xr2 = x2 * r2i; yr2 = y2 * r2i; sour = -dl2*(log(r2) + d12 * dr2 *(yr2 -xr2)); dipl = y*dl*r2i*(1. +d12*dr2*(4.0*xr2 - 1.0)); else % Near Field Closed Form Expressions ier = 0; xx1 = x + dl2; xx2 = x - dl2; r1 = sqrt(xx1*xx1 + y2); r2 = sqrt(xx2*xx2 + y2); if (y2 > 0.25e-12); dipl = -sign(y)*(acos(xx1/r1) - acos(xx2/r2)); end axx1 = abs(xx1); axx2 = abs(xx2); if (axx1 < tol); dl1 = 0.0; else dl1 = xx1*log(r1); end if (axx2 < tol); dl2 = 0.0; else dl2 = -xx2*log(r2); end sour = -(y*dipl +dl1 + dl2 -dl); end