三维结构化网格里面的六面体单元任一个“面”上的四个点一定共面吗？

flyboys · 发表于 2005-6-25 00:55:06

马上注册，结交更多好友，享用更多功能，让你轻松玩转社区。

您需要登录才可以下载或查看，没有账号？注册

x

[这个贴子最后由flyboys在 2005/06/25 12:46pm 第 1 次编辑]

三维结构化网格里面的六面体单元任一个“面”上的四个点一定共面吗？
如果不共面，一般是如何得到这个面的面积矢量，即面积和其法矢量的乘积？
弄不明白，请有个编程经验的大侠告诉我！谢谢

harry · 发表于 2005-6-25 18:52:10

   INCLUDE ';COM3D1.FI';
C
      MFLAG1=0
      MFLAG2=0
      MFLAG3=0
      MFLAG4=0
      MFLAG5=0
      MFLAG6=0
   IF(    IFLAG.EQ.1) THEN
      MFLAG1=1
   ELSE IF(IFLAG.EQ.2) THEN
      MFLAG2=1
   ELSE IF(IFLAG.EQ.3) THEN
      MFLAG3=1
   ELSE IF(IFLAG.EQ.4) THEN
      MFLAG4=1
   ELSE IF(IFLAG.EQ.5) THEN
      MFLAG5=1
   ELSE IF(IFLAG.EQ.6) THEN
      MFLAG6=1
   ENDIF
C
      MWAL = MWALN*MFLAG1 + MWALS*MFLAG2 + MWALE*MFLAG3
   >    + MWALW*MFLAG4 + MWALT*MFLAG5 + MWALB*MFLAG6
   ISANOW =  ISAN*MFLAG1 +  ISAS*MFLAG2 +  ISAE*MFLAG3
   >    +  ISAW*MFLAG4 +  ISAT*MFLAG5 +  ISAB*MFLAG6
   DO 1110 IM=1,MWAL
      IIP= INDWN(IM)*MFLAG1 + INDWS(IM)*MFLAG2 + INDWE(IM)*MFLAG3
   >    + INDWW(IM)*MFLAG4 + INDWT(IM)*MFLAG5 + INDWB(IM)*MFLAG6
   IIB=IIP+       NI*MFLAG1 -       NI*MFLAG2 +          MFLAG3
   >    -          MFLAG4 +    NIJ*MFLAG5 -    NIJ*MFLAG6
   ID1=IIP+(1+NI+NIJ)*MFLAG1 +    NIJ*MFLAG2 + (1+NIJ)*MFLAG3
   >    +  (NI+NIJ)*MFLAG4 +  (NI+NIJ)*MFLAG5 + (NI+1)*MFLAG6
   ID2=IIP+       NI*MFLAG1 +          MFLAG2 + (NI+1)*MFLAG3
   >                         + (NIJ+1)*MFLAG5
   ID3=IIP+  (NI+NIJ)*MFLAG1 + (NIJ+1)*MFLAG2 +(1+NI+NIJ)*MFLAG3
   >    +    NIJ*MFLAG4 +(1+NI+NIJ)*MFLAG5 +       NI*MFLAG6
   ID4=IIP+ (NI+1)*MFLAG1                   +          MFLAG3
   >    +       NI*MFLAG4 +    NIJ*MFLAG5 +          MFLAG6
C    NORMAL VECTOR, AREA AND NORMAL DISTANCE
         DAX=F(ISCOX+ID1)-F(ISCOX+ID2)
         DAY=F(ISCOY+ID1)-F(ISCOY+ID2)
         DAZ=F(ISCOZ+ID1)-F(ISCOZ+ID2)
         DBX=F(ISCOX+ID3)-F(ISCOX+ID4)
         DBY=F(ISCOY+ID3)-F(ISCOY+ID4)
         DBZ=F(ISCOZ+ID3)-F(ISCOZ+ID4)
      DN1=DAY*DBZ-DBY*DAZ
      DN2=DBX*DAZ-DAX*DBZ
      DN3=DAX*DBY-DAY*DBX
      AREA=0.5*SQRT(DN1*DN1+DN2*DN2+DN3*DN3)
      A1=0.5*DN1/(AREA+SMALL)
      A2=0.5*DN2/(AREA+SMALL)
      A3=0.5*DN3/(AREA+SMALL)
   DELTA=ABS((F(ISX+IIP)-F(ISX+IIB))*A1+(F(ISY+IIP)-F(ISY+IIB))*A2
   >       +(F(ISZ+IIP)-F(ISZ+IIB))*A3)+SMALL
C    SHEAR STRESS AT THE WALL
         A11=1.-A1*A1
         A22=1.-A2*A2
         A33=1.-A3*A3
         A12=A1*A2
         A13=A1*A3
         A23=A2*A3
      VXCOMP=A11*F(ISU+IIP)-A12*F(ISV+IIP)-A13*F(ISW+IIP)
      VYCOMP=A22*F(ISV+IIP)-A12*F(ISU+IIP)-A23*F(ISW+IIP)
      VZCOMP=A33*F(ISW+IIP)-A13*F(ISU+IIP)-A23*F(ISV+IIP)
      VPARL=SQRT(VXCOMP*VXCOMP+VYCOMP*VYCOMP+VZCOMP*VZCOMP)
      SQRTK=SQRT(ABS(F(ISTE+IIP)))
      YPLUS=F(ISDEN+IIP)*SQRTK*CDQR*DELTA/VISCOS
      TMULT=VISCOS/DELTA
      IF(YPLUS.GT.11.6)
   > TMULT=F(ISDEN+IIP)*CDQR*SQRTK*CAPPA/ALOG(ELOG*YPLUS)
C    MODIFICATION OF SOURCE TERMS
      FACTOR=TMULT*AREA
      UPUB=F(ISU+IIP)-F(ISU+IIB)
      VPVB=F(ISV+IIP)-F(ISV+IIB)
      WPWB=F(ISW+IIP)-F(ISW+IIB)
   IF( NN.EQ.2) THEN
   F(ISSP+IIP)=F(ISSP+IIP)-FACTOR*A11
   F(ISSU+IIP)=F(ISSU+IIP)+FACTOR*(F(ISU+IIB)*A11+VPVB*A12+WPWB*A13)
   ELSEIF(NN.EQ.3) THEN
   F(ISSP+IIP)=F(ISSP+IIP)-FACTOR*A22
   F(ISSU+IIP)=F(ISSU+IIP)+FACTOR*(F(ISV+IIB)*A22+UPUB*A12+WPWB*A23)
   ELSEIF(NN.EQ.4) THEN
   F(ISSP+IIP)=F(ISSP+IIP)-FACTOR*A33
   F(ISSU+IIP)=F(ISSU+IIP)+FACTOR*(F(ISW+IIB)*A33+UPUB*A13+VPVB*A23)
   ELSEIF(NN.EQ.5) THEN
      IF(YPLUS.LE.11.6) THEN
         TERM=F(ISDEN+IIP)*CDTQ*SQRTK*YPLUS/DELTA
         ELSE
         TERM=F(ISDEN+IIP)*CDTQ*SQRTK*ALOG(ELOG*YPLUS)/(CAPPA*DELTA)
      ENDIF
      TAUW=TMULT*VPARL
      GENR=TAUW*TAUW/(CAPPA*VISCOS*YPLUS+SMALL)
   F(ISSP+IIP)=TESPC(IIP)-TERM*F(ISVOLP+IIP)
   F(ISSU+IIP)=TESUC(IIP)+GENR*F(ISVOLP+IIP)
   ELSEIF(NN.EQ.6) THEN
   F(ISSP+IIP)=-GREAT
   F(ISSU+IIP)= GREAT*ABS(F(ISTE+IIP))**1.5*CDTQ/(CAPPA*DELTA)
CTANG.TANG.TANG.--------------------------------------------------
   ELSEIF(NN.EQ.7) THEN
   F(ISS+IIB)=F(ISS+IIP)
CTANG.TANG.TANG.--------------------------------------------------
   ENDIF
1110 F(ISANOW+IIP)=0.0
   RETURN
   END

flyboys · 发表于 2005-6-26 14:53:51

老兄，没有程序说明看不懂啊，好象还是有限元程序！

anyone · 发表于 2005-6-28 00:36:51

四个点不一定共面. 算法可以有不同的, 比如把面分成3角形. 重点是唯一性(每个面总用同一个算法)和一至性(所以面都用类似算法).
一般做法: 先算面积矢量, 也就是面在3个方向投影, 这是不变量, 和3角形划分无关. 面积可以用面积矢量的模.

		自动登录	找回密码
密码			注册

三维结构化网格里面的六面体单元任一个“面”上的四个点一定共面吗？

马上注册，结交更多好友，享用更多功能，让你轻松玩转社区。