CFD Online Logo CFD Online URL
www.cfd-online.com
[Sponsors]
Home > Forums > Main CFD Forum

Can anyone explain about this function

Register Blogs Members List Search Today's Posts Mark Forums Read

Reply
 
LinkBack Thread Tools Display Modes
Old   October 3, 2010, 04:41
Default Can anyone explain about this function
  #1
New Member
 
cxcxcx0505
Join Date: May 2010
Posts: 17
Rep Power: 7
cxcxcx0505 is on a distinguished road
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
Attached Images
File Type: jpg untitled.jpg (100.9 KB, 6 views)
cxcxcx0505 is offline   Reply With Quote

Reply

Thread Tools
Display Modes

Posting Rules
You may not post new threads
You may not post replies
You may not post attachments
You may not edit your posts

BB code is On
Smilies are On
[IMG] code is On
HTML code is Off
Trackbacks are On
Pingbacks are On
Refbacks are On


Similar Threads
Thread Thread Starter Forum Replies Last Post
Compile problem ivanyao OpenFOAM Running, Solving & CFD 1 October 12, 2012 09:31
BlockMesh FOAM warning gaottino OpenFOAM Native Meshers: blockMesh 7 July 19, 2010 14:11
latest OpenFOAM-1.6.x from git failed to compile phsieh2005 OpenFOAM Bugs 25 February 9, 2010 05:37
Error with Wmake skabilan OpenFOAM Installation 3 July 28, 2009 00:35
Axisymmetrical mesh Rasmus Gjesing (Gjesing) OpenFOAM Native Meshers: blockMesh 10 April 2, 2007 14:00


All times are GMT -4. The time now is 00:44.