CFD Online Discussion Forums

CFD Online Discussion Forums (https://www.cfd-online.com/Forums/)
-   Main CFD Forum (https://www.cfd-online.com/Forums/main/)
-   -   Can anyone explain about this function (https://www.cfd-online.com/Forums/main/80674-can-anyone-explain-about-function.html)

cxcxcx0505 October 3, 2010 04:41

Can anyone explain about this function
 
1 Attachment(s)
Do anyone know about this metric function?Where can I find more information about this?And,I have another question,what is the fmag and fangle for in this function?
Thank you.


Subroutine Metric

Use young

Integer :: index(4)
Real :: dum1,dum2,dum3,dum4
Real,Dimension(mx,my) :: xxi,xeta,yxi,yeta,xtau,ytau
Real :: fudh,fvdh,fanglet
Real :: fxtemp1,fytemp1,fxtemp2,fytemp2
Real :: xvel,yvel,fcfind

index(1)=ios-1; index(2)=ioe; index(3)=jos-1; index(4)=joe

Do j=jos,joe-1
Do i=ios-1,ioe
xtau(i,j)=(nxc(i,j)-xc(i,j))/dt
ytau(i,j)=(nyc(i,j)-yc(i,j))/dt
End Do
End Do

Call Compact2D_xi(xc,xxi,mx,my,index)
Call Compact2D_xi(yc,yxi,mx,my,index)
Call Compact2D_eta(xc,xeta,mx,my,index)
Call Compact2D_eta(yc,yeta,mx,my,index)

Do j=jos-1,joe
Do i=ios-1,ioe
vol(i,j)=xxi(i,j)*yeta(i,j)-xeta(i,j)*yxi(i,j)
xit(i,j)=(-xtau(i,j)*yeta(i,j)+xeta(i,j)*ytau(i,j))/vol(i,j)
etat(i,j)=(xtau(i,j)*yxi(i,j)-xxi(i,j)*ytau(i,j))/vol(i,j)
xix(i,j)=yeta(i,j)/vol(i,j)
xiy(i,j)=-xeta(i,j)/vol(i,j)
etax(i,j)=-yxi(i,j)/vol(i,j)
etay(i,j)=xxi(i,j)/vol(i,j)
q11(i,j)=(xix(i,j)**2+xiy(i,j)**2)*vol(i,j)
q22(i,j)=(etax(i,j)**2+etay(i,j)**2)*vol(i,j)
End Do
End Do

Do j=js,je
Do i=is,ie
x1=xc(i-1,j-1); x2=xc(i,j-1); x3=xc(i,j); x4=xc(i-1,j)
y1=yc(i-1,j-1); y2=yc(i,j-1); y3=yc(i,j); y4=yc(i-1,j)
xx=x(i,j); yy=y(i,j)

dx1=xx-x1; dx2=xx-x2; dx3=xx-x3; dx4=xx-x4
dy1=yy-y1; dy2=yy-y2; dy3=yy-y3; dy4=yy-y4

xi1=xix(i-1,j-1)*dx1+xiy(i-1,j-1)*dy1
eta1=etax(i-1,j-1)*dx1+etay(i-1,j-1)*dy1

xi2=1.+xix(i,j-1)*dx2+xiy(i,j-1)*dy2
eta2=etax(i,j-1)*dx2+etay(i,j-1)*dy2

xi3=1.+xix(i,j)*dx3+xiy(i,j)*dy3
eta3=1.+etax(i,j)*dx3+etay(i,j)*dy3

xi4=xix(i-1,j)*dx4+xiy(i-1,j)*dy4
eta4=1.+etax(i-1,j)*dx4+etay(i-1,j)*dy4

xic=0.25*(xi1+xi2+xi3+xi4)
etac=0.25*(eta1+eta2+eta3+eta4)

rr(i,j,1)=(1.-etac)*(1.-xic)
rr(i,j,2)=(1.-etac)*xic
rr(i,j,3)=etac*(1.-xic)
rr(i,j,4)=etac*xic
End Do
End Do

! Compute x- and y- velocities of each point
fudh=(fmag2-fmag1)/dt
fvdh=0.
fanglet=(fangle2-fangle1)/dt

Do j=jos-1,joe
Do i=ios-1,ioe
fxtemp1=xc(i,j)*cos(stroke)+yc(i,j)*sin(stroke)-fmag1
fytemp1=-xc(i,j)*sin(stroke)+yc(i,j)*cos(stroke)
fxtemp2=fxtemp1*cos(fangle1)+fytemp1*sin(fangle1)
fytemp2=-fxtemp1*sin(fangle1)+fytemp1*cos(fangle1)
fcfind=fxtemp2**2/(aa**2)+fytemp2**2/(bb**2)

If (fcfind <= 1.) Then
xvel=fudh-fxtemp1*fanglet*sin(fangle2-fangle1)-fytemp1*fanglet*cos(fangle2-fangle1)
yvel=fvdh+fxtemp1*fanglet*cos(fangle2-fangle1)-fytemp1*fanglet*sin(fangle2-fangle1)
uob(i,j)=xvel*cos(stroke)-yvel*sin(stroke)
vob(i,j)=xvel*sin(stroke)+yvel*cos(stroke)
ck(i,j)=ckzero
Else
uob(i,j)=0.
vob(i,j)=0.
ck(i,j)=ckinf
End If
End Do
End Do

cmapxi=0.
cmapeta=0.

Do j=jos,joe-1
Do i=ios,ioe-1
xgrad=abs(ck(i+1,j)-2.*ck(i,j)+ck(i-1,j))
ygrad=abs(ck(i,j+1)-2.*ck(i,j)+ck(i,j-1))
If (xgrad > 1.) Then
cmapxi(i+1,j)=1.; cmapxi(i,j)=1.; cmapxi(i-1,j)=1.
End If
If (ygrad > 1.) Then
cmapeta(i,j+1)=1.; cmapeta(i,j)=1.; cmapeta(i,j-1)=1.
End If
End Do
End Do

Return
End Subroutine Metric


All times are GMT -4. The time now is 12:42.