Navigation

    CFD中文网

    CFD中文网

    • Login
    • Search
    • 最新
    1. Home
    2. 程迪
    程
    • Profile
    • Following
    • Followers
    • Topics
    • Posts
    • Groups

    猫狗大爷

    @程迪

    航天十一院工程师,UConn访学ing

    186
    Posts
    2298
    Profile views
    20
    Followers
    0
    Following
    Joined Last Online
    Website chengdi123000.github.io Location University of Connecticut Age 32

    程迪 Follow

    Best posts made by 程迪

    • RE: openfoam中cyclic周期性边界的问题

      cyclic BC是个大坑,因为OF的理论和编程经常不自洽。感觉网上从来没说清楚过。刚才我才搞明白。

      • 首先,face addressing的owner, neighbour和ldu addressing的lower, upper并不天然等价。face addressing描述的的网格的几何拓扑关系。而ldu addressing需要描述的的是数值拓扑关系,如果没有耦合边界,二者是等价的,但是如果有耦合边界,就有了face addressing描述不了的off-diagonal term。这就有两种处理方式:

      • 在OF+1706中,有两种lduAddressing,包括fvMeshLduAddressing和fvMeshPrimitiveLduAddressing。

        • 实际上fvMeshLduAddressing是将face addressing和ldu addressing等价了,内部面的owner就是lower,neighbour就是upper。这也是最最常用的类型,可以称之为普通ldu。
        • 而fvMeshPrimitiveLduAddressing(这里的Primitive的意思是not developed or derived from anything else )中的lower可以添加owner之外的相互作用关系。但是这个类位于overset目录下,可见它主要用于重叠网格的情况,可以称之为overset ldu。
      • 而coupled BC是实现在普通ldu中。依靠两个类:coupledFvPatch和coupledFvPatchField,其中关键是他们各自继承了lduInterface和lduInterfaceField,并实现了相关函数。所以其实interface和coupled BC是几乎同样的概念。

      • 从概念上讲,对于普通LDU,在有interface的情况下虽然在lduMatrix中存了upper,diag,lower矩阵系数,但是还有额外耦合项,也就是说$A=L+D+U+C\ne L+D+U$,$C$就是额外耦合项,由于数据结构设计上的限制,不能存到普通ldu的矩阵中,需要单独处理。

      • 处理方式也比较简单,$C$因为也是代表cell-cell相互作用,所以主要是非对角的,因为对角的部分完全可以包含在$D$中,而且源项部分也可以容纳在原来的源项里,可以分裂为上下两部分,$C=C_L+C_U$。在lduMatrix.Amul()的操作中,计算$A\cdot x$的操作被分为$D\cdot x+(L+U)\cdot x +C_U\cdot x_n $:

        • initMatrixInterfaces: $C_L\cdot x_o$,其中$x_o$是$x$中属于耦合边界的内侧单元,计算结果要送到另一侧去;
        • all cells: ApsiPtr[cell] = diagPtr[cell]*psiPtr[cell];: $D\cdot x$
        • all faces: ApsiPtr[uPtr[face]] += lowerPtr[face]*psiPtr[lPtr[face]];: $L\cdot x$
        • all faces: ApsiPtr[lPtr[face]] += upperPtr[face]*psiPtr[uPtr[face]];: $U\cdot x$
        • updateMatrixInterfaces: $C_U\cdot x_n$,其中$x_n$只包含$x$中属于耦合边界的外侧单元;
        • 这种$C\cdot x$叫做patch Contribution。 $C_L, C_U$叫interfaceUpper, interfaceLower。
      • 问:为什么$C_L,C_U$要分开处理?答:因为在存在processor边界时,$x_n$并不在本地,要通过通信获取结果,或者获取系数。

      • 问:$C_L$和$C_U$存在哪儿?答:

        • class fvMatrix
          :
              public refCount,
              public lduMatrix
          {
          //...
                  //- Boundary scalar field containing pseudo-matrix coeffs
                  //  for internal cells
                  FieldField<Field, Type> internalCoeffs_; //注意这里不是引用,是真货
          
                  //- Boundary scalar field containing pseudo-matrix coeffs
                  //  for boundary cells
                  FieldField<Field, Type> boundaryCoeffs_; //都初始化为0;
          //...
          }
          
          
      • coupled BC, interface, constraint BC等的关系:
        interface就是coupled BC,constraint BC是OF文档中的分类,其实包含了empty, wedge,symmetry等几何约束类的BC和coupled BC。coupled BC中比较重要的两类是cyclic和processor。

      总的来看cyclic边界条件本身的实现应该没有问题,还挺巧妙的。

      回到这个问题,我觉得可能是CFL数计算的部分可能在cyclic BC处有bug。
      Courant数计算方式来看,普通Courant数和其他一样,额外多了一个Interface Courant数。里面多了个这玩意儿:

          scalarField sumPhi
          (
              mixture.nearInterface()().primitiveField()
             *fvc::surfaceSum(mag(phi))().primitiveField()
          );
      
          alphaCoNum = 0.5*gMax(sumPhi/mesh.V().field())*runTime.deltaTValue();
      
          meanAlphaCoNum =
              0.5*(gSum(sumPhi)/gSum(mesh.V().field()))*runTime.deltaTValue();
      
      // mixture.nearInterface()(),定义是:
      
       Foam::tmp<Foam::volScalarField>
      Foam::interfaceProperties::nearInterface() const
      {
           return pos(alpha1_ - 0.01)*pos(0.99 - alpha1_);
      }
      

      可能是相边界速度太大了吧。估算一下有上万了,你的平均Courant这么小,最大Courant这么大。。。你又是这么均匀的网格。。

      posted in OpenFOAM
      程
      程迪
    • RE: 进口马赫数0.7左右,用rhosimpleFoam不收敛的问题

      @小葱小虫

      要不算个初始场?

      在fvSolution中加

      solvers
      {
      ...
           Phi
           {
                solver PCG;
                preconditioner PIC;
                tolerance 1e-12;
                relTol 0;
           }
      }
      
      potentialFlow
      {
          nCorrectors 10;
      }
      

      然后运行potentialFlow搞个速度初始场试试?

      或者在controlDict中加一些debugSwitches

      DebugSwitches
      {
      lduMatrix 2;
      ...// 去$FOAM_ETC/controlDict下找和湍流有关的。
      }
      

      还有就是加functionObject : fieldMinMax看看是哪个位置的场出了问题

          fieldMinMax1
          {
              type        fieldMinMax;
              libs ("libfieldFunctionObjects.so");
              write       yes;
              log         yes;
              location    yes; // 关键是这个!
              mode        magnitude;
              fields
              (
                  U
                  p
                  omega
              );
          }
      
      posted in OpenFOAM
      程
      程迪
    • RE: solveSegragated的问题?

      哦,

       diag() = saveDiag;//恢复diag
      

      是在括号之前的,所以每次的修改不会叠加。

      之前理解错误的一点是所有的边界有关的系数都是放在internalCoeffs和boundaryCoeffs中的,不仅仅是coupled BC。从初始化时候的size可以看出来。

      template<class Type>
      Foam::fvMatrix<Type>::fvMatrix
      (
          const GeometricField<Type, fvPatchField, volMesh>& psi,
          const dimensionSet& ds
      )
      :
          lduMatrix(psi.mesh()),
          psi_(psi),
          dimensions_(ds),
          source_(psi.size(), Zero),
          internalCoeffs_(psi.mesh().boundary().size()), //和边界的face数大小一样。
          boundaryCoeffs_(psi.mesh().boundary().size()),
          faceFluxCorrectionPtr_(nullptr)
      {
      //...
      

      不过非耦合的应该是internalCoeffs放比例系数,boundaryCoeffs放源项,而耦合的BC是internalCoeffs放owner的系数,boundaryCoeffs放的另一侧单元的系数。从addBoundarySource()的源代码可以看出。

      
      template<class Type>
      void Foam::fvMatrix<Type>::addBoundarySource
      (
          Field<Type>& source,
          const bool couples
      ) const
      {
          forAll(psi_.boundaryField(), patchi)
          {
              const fvPatchField<Type>& ptf = psi_.boundaryField()[patchi];
              const Field<Type>& pbc = boundaryCoeffs_[patchi];
      
              if (!ptf.coupled())//非耦合边界,只有owner,另一侧没有单元,pbc存的是边界源项。
              {
                  addToInternalField(lduAddr().patchAddr(patchi), pbc, source);
              }
              else if (couples)//耦合边界,如果couples==true, 另一侧有单元,pbc存的是另一侧单元的系数。
              {
                  const tmp<Field<Type>> tpnf = ptf.patchNeighbourField();
                  const Field<Type>& pnf = tpnf();
      
                  const labelUList& addr = lduAddr().patchAddr(patchi);
      
                  forAll(addr, facei)
                  {
                      source[addr[facei]] += cmptMultiply(pbc[facei], pnf[facei]);
                  }
              }
          }
      }
      
      posted in OpenFOAM
      程
      程迪
    • RE: 有关low Mach的疑问

      @李东岳

      其实密度基/压力基,耦合/分离,显式/隐式 是几个正交的形容词。

      比如foam_extend的coupledPU就是耦合的压力基
      常规OF里的pimpleFoam就是分离的压力基

      DensityBasedTurbo就是显式的耦合的密度基
      OF的rhoCentralFoam就是显式的分离的密度基

      密度基和压力基求解压力时用的压力泊松方程/压力修正方程 还是 状态方程,压力泊松/修正方程是微分方程,状态方程式代数方程,所以差别比较大。

      耦合/分离看你的主变量在每个迭代/时间步里是同时更新还是分别更新的。

      隐式/显式看你的空间离散用没用未知时间步的量。具体隐式还分点隐/全隐之类的玩意儿。甚至可以半隐式。

      posted in OpenFOAM
      程
      程迪
    • RE: 有关low Mach的疑问

      没有theory guide的弊端。。。

      posted in OpenFOAM
      程
      程迪
    • RE: OpenFOAM不同版本代码转换的问题

      据我经验,错误最多的是这么几个:

      1. of4 的const ref 和ref分得比较清楚了,参考http://www.openfoam.com/documentation/developer-upgrade-guide.php
      2. of版本升级时,热力学那块类变化比较大
      3. of中cxx类有abstract base class (ABC), ABC不能实例化,需要在子类中实现所有的虚方法,但是of版本升级时可能虚方法数量发生了变化,这个在错误中会有提示。
      4. Make/options里include path变化,特别湍流/热力学那些库变化比较大
      纯虚函数个数 of41 of2x
      fvcDdt 5 4
      fvmDdt 4 3
      fvcDdtPhiCorr 2 2
      fvcDdtUfCorr 2 0
      meshPhi 1 1
      type 1 1

      不过of41支持cxx11了,似乎有些语法糖衣可以用:比如auto x = phi.ref(),不用写又臭又长的模板展开了。

      posted in OpenFOAM
      程
      程迪
    • RE: CFL3D开源了

      安装脚本
      大家可以试着玩一玩

      # MPI和fortran,UBUNTU为例
      sudo apt-get update
      sudo apt-get install libopenmpi-dev gfortran build-essential git -y
      # 下载chengdi123000修正过的cfl3d代码,(spalart.F, cputim.F)
      cd $HOME
      mkdir cfl3d
      cd cfl3d
      git clone https://github.com/chengdi123000/CFL3D.git
      # cgns代码和安装
      cd $HOME/cfl3d
      wget "http://jaist.dl.sourceforge.net/project/cgns/cgnslib_2.5/Release 5/cgnslib_2.5-5.tar.gz" -O cgnslib_2.5-5.tar.gz
      tar xf cgnslib_2.5-5.tar.gz
      mkdir cgns
      mkdir cgns/include
      mkdir cgns/lib
      cd cgnslib_2.5-5
      ./configure --prefix=../cgns #这是按照最简单的配置进行编译的,没有配hdf5,并行之类高级的东西。
      make
      make install
      # 生成makefile
      cd $HOME/cfl3d/CFL3D/build
      ./Install -fastio -cgnsdir=$HOME/cfl3d/cgns
      # 编译cfl3d
      cd $HOME/cfl3d/CFL3D/build
      cp makefile_linux_gfortran_openmpi makefile #覆盖掉新生成的makefile
      make cfl3d_seq cfl3d_mpi splitter cfl3d_tools
      make precfl3d
      make ronnie preronnie
      make maggie
      make splittercmplx cfl3dcmplx_seq cfl3dcmplx_mpi 
      # 建立链接
      cd $HOME/cfl3d/CFL3D
      mkdir bin
      cd bin
      ## basic
      ln -s ../build/cfl/seq/cfl3d_seq
      ## parallel version
      ln -s ../build/cfl/mpi/cfl3d_mpi
      ## mesh block splitter
      ln -s ../build/split/seq/splitter
      ## memory usage estimator
      ln -s ../build/precfl/seq/precfl3d
      ## overset mesh related tool
      ln -s ../build/mag/seq/maggie
      ## mesh deformation tool and its memory usage estimator
      ln -s ../build/ron/seq/ronnie
      ln -s ../build/preron/seq/preronnie
      ## complex version of cfl3d which is used to compulate flight derivative
      ln -s ../build/cflcmplx/seq/cfl3dcmplx_seq
      ln -s ../build/cflcmplx/mpi/cfl3dcmplx_mpi
      ln -s ../build/splitcmplx/seq/splittercmplx
      ## cfl3d_tools
      ln -s ../build/tools/seq/v6inpdoubhalf
      ln -s ../build/tools/seq/Get_FD
      ln -s ../build/tools/seq/initialize_field
      ln -s ../build/tools/seq/nmf_to_cfl3dinput
      ln -s ../build/tools/seq/cgns_to_cfl3dinput
      ln -s ../build/tools/seq/v6inpswitchijk
      ln -s ../build/tools/seq/v6_ronnie_mod
      ln -s ../build/tools/seq/moovmaker
      ln -s ../build/tools/seq/grid_perturb_cmplx
      ln -s ../build/tools/seq/v6_restart_mod
      ln -s ../build/tools/seq/gridswitchijk
      ln -s ../build/tools/seq/cfl3d_to_nmf
      ln -s ../build/tools/seq/cgns_readhist
      ln -s ../build/tools/seq/p3d_to_cfl3drst
      ln -s ../build/tools/seq/everyother_xyz
      ln -s ../build/tools/seq/p3d_to_INGRID
      ln -s ../build/tools/seq/v6_ronnie_mod.F90
      ln -s ../build/tools/seq/INGRID_to_p3d
      ln -s ../build/tools/seq/grid_perturb
      ln -s ../build/tools/seq/cfl3dinp_to_FVBND
      ln -s ../build/tools/seq/cfl3d_to_pegbc
      ln -s ../build/tools/seq/plot3dg_to_cgns
      ln -s ../build/tools/seq/XINTOUT_to_ovrlp
      # 下载算例
      cd $HOME/cfl3d
      mkdir TestCases
      cd TestCases
      wget https://cfl3d.larc.nasa.gov/Cfl3dv6/2DTestcases/Flatplate/Flatplate.tar.Z
      wget https://cfl3d.larc.nasa.gov/Cfl3dv6/2DTestcases/Flatplateskew/Flatplateskew.tar.Z
      wget https://cfl3d.larc.nasa.gov/Cfl3dv6/2DTestcases/Flatplateyplus/Flatplateyplus.tar.gz
      wget https://cfl3d.larc.nasa.gov/Cfl3dv6/2DTestcases/Backstep/Backstep.tar.gz
      wget https://cfl3d.larc.nasa.gov/Cfl3dv6/2DTestcases/Transdiff/Transdiff.tar.Z
      wget https://cfl3d.larc.nasa.gov/Cfl3dv6/2DTestcases/NACA_4412/NACA_4412.tar.Z
      wget https://cfl3d.larc.nasa.gov/Cfl3dv6/2DTestcases/RAE_Sensitivity/RAE_Sensitivity.tar.gz
      wget https://cfl3d.larc.nasa.gov/Cfl3dv6/2DTestcases/Ramp/Ramp.tar.Z
      wget https://cfl3d.larc.nasa.gov/Cfl3dv6/2DTestcases/Cylinder/Timeaccstudy.tar.Z
      wget https://cfl3d.larc.nasa.gov/Cfl3dv6/2DTestcases/N0012/Spaceaccstudy.tar.Z
      wget https://cfl3d.larc.nasa.gov/Cfl3dv6/2DTestcases/Ejectornozzle/Ejectornozzle.tar.gz
      wget https://cfl3d.larc.nasa.gov/Cfl3dv6/2DTestcases/Pitch/Pitch0012.tar.Z
      wget https://cfl3d.larc.nasa.gov/Cfl3dv6/2DTestcases/Rotorstator/Rotorstator.tar.Z
      wget https://cfl3d.larc.nasa.gov/Cfl3dv6/2DTestcases/Hump/Humpcase.tar.gz
      wget https://cfl3d.larc.nasa.gov/Cfl3dv6/2DTestcases/Curvature/SoMellor.tar.gz
      wget https://cfl3d.larc.nasa.gov/Cfl3dv6/3DTestcases/Axibump/Axibump.tar.gz
      wget https://cfl3d.larc.nasa.gov/Cfl3dv6/3DTestcases/ONERA_M6/ONERA_M6.tar.Z
      wget https://cfl3d.larc.nasa.gov/Cfl3dv6/3DTestcases/ARA_M100/ARA_M100.tar.Z
      wget https://cfl3d.larc.nasa.gov/Cfl3dv6/3DTestcases/ARA_M100_XMERA/ARA_M100_XMERA.tar.Z
      wget https://cfl3d.larc.nasa.gov/Cfl3dv6/3DTestcases/Delta/Delta_cgns.tar.Z
      
      posted in OpenFOAM
      程
      程迪
    • RE: AeroFoamV1-2008 移植到OF41

      @李东岳

      说得跟白马非马一样。。。

      FVM里我的理解是重构=插值,FVM等价于二阶Modal DG

      rCF 的KT/KNP格式其实也就是进化evolution步,总之是要算出数值通量的那一步。

      而且rCF的其实不是教科书上那种密度基,可以算是分离的,segregated的密度基,非常interesting。

      你看rCF的代码,丫的先算出密度,在算动量,然后用动量除以密度算速度,最后算总能,用总能除以密度减去动能得比内能,再和密度一起算压力。

      而AeroFoam,DBT求解器都是先算出守恒量再统一算出原始变量的。

      其实二者在密度和动量方程上差不多,但是在能量方程那一项和粘性项差别就很大了。AF和DBT的代码n+1时刻的变量只直接依赖于n时刻的变量。但是rCF的n+1时刻的密度和速度只依赖于n时刻的值,n+1时刻的能量和粘性项直接依赖的是n+1时刻的速度和密度。

      这种搞法不知道对不对,但是似乎大概信息传递能快一点点。代价是不一致性,也就是rCF在n+1时刻的守恒量和原始变量是不对应的。

      posted in OpenFOAM
      程
      程迪
    • RE: OF可压流求解器

      https://github.com/chengdi123000/compressibleFoam/tree/of41
      把 http://pavanakumar.github.io/compressibleFoam 移植到了openfoam 4.1下。里面有roe格式。没啥大区别。

      posted in OpenFOAM
      程
      程迪
    • RE: sutherland和janaf

      @CC_Cherry

      https://www.cfd-online.com/Wiki/Sutherland's_law

      但是Sutherland's Law 在https://turbmodels.larc.nasa.gov/implementrans.html上还有另一种形式,两者是可以相互转化的,参考

      https://en.wikipedia.org/wiki/Viscosity#Effect_of_temperature_on_the_viscosity_of_a_gas

      不过,NASA TMR给的空气粘性和网上各个渠道给的有些不一样。我算出来的是这样的。

      输运性质

      根据NASA TMR网站的说明,空气粘性为:
      \begin{equation}
      \mu = \mu_0 ( \frac{T}{T_0} )^{3/2} ( \frac{T_0+S}{T+S} )$
      \end{equation}
      其中,

      • $\mu_0$ = 1.716E-5 kg/m/s;
      • $T_0$ = 491.6 degR = 273.111111 K;
      • $S$= 198.6 degR = 110.333333 K;

      并要求:Pr = 0.72, 湍流Prt = 0.9

      但是OpenFOAM中的sutherland模型是:
      \begin{equation}
      \mu = \frac{A_s\sqrt T}{1+T_s/T}
      \end{equation}
      因此,从NASA给定的粘性可以计算出:

      • $A_s = \frac{\mu_0(T_0+S)}{T_0^{3/2}}$ = 1.4578427435006915e-06 kg/m/s/sqrt(K);
      • Ts = S = 110.333333 K;

      从而可得系数:

      transport //sutherland model
      {
      As 1.4578427435006915e-06;
      Ts 110.333333;
      }
      

      但是,OpenFOAM中sutherland模型不能指定Pr数,而是通过如下代码计算的:

      // @ src/thermophysicalModels/specie/transport/sutherland/sutherlandTransportI.H
      
      template<class Thermo>
      inline Foam::scalar Foam::sutherlandTransport<Thermo>::kappa
      (
          const scalar p, const scalar T
      ) const
      {
          scalar Cv_ = this->Cv(p, T);
          return mu(p, T)*Cv_*(1.32 + 1.77*this->R()/Cv_);
      }
      

      这等价于
      \begin{equation}
      Pr = \frac{\gamma}{1.32+1.77(\gamma-1)}
      \end{equation}
      对于gamma=1.4的情形,Pr=0.6903353057199211,因此会高估热导率

      根据OpenFOAM的bug report,这个问题会在OpenFOAM 4.x和OpenFOAM dev中得到修复。

      posted in OpenFOAM
      程
      程迪

    Latest posts made by 程迪

    • RE: 最快速的上google的方式是什么?

      V2ray自建服务器,有一键脚本,联通不封github v2ray的repo,可以直接下载。手机用v2rayng,电脑用v2rayn。

      服务器可选在日本、荷兰或者美国,最好根据当地联通/电信/移动的网络来选。

      反正,u2b 4k能看。而且还能做透明代理啥的。也可以做订阅。

      posted in C斯达克
      程
      程迪
    • RE: 压力方程和解

      这个是理论大坑啊。

      如果压力基的求解器,求解的是速度和压力场,那么得用连续性方程来求压力。
      如果是密度基的,连续性可以求密度,然后用状态方程求压力。

      posted in Algorithm
      程
      程迪
    • RE: CPU在并行时候的指派

      Linux调度的话默认多核是会给任务分配不同的核直到所有核都满。但是你调整了优先级,那个NI,nice value。

      posted in Algorithm
      程
      程迪
    • RE: 求助如何在HPC Cluster上安装OF5.0
      1. 最好是让管理员装。
      2. 如果比较新的Linux,有singular或者类似的docker容器技术,可以直接把镜像弄上去就能用。
      3. 如果没有,可以用HPC包管理器,自动安装和解决依赖,比如spack 或者easy build。一个是美国货,一个欧洲货,建议采用美国货。如果能联网,这是最快最方便的。
      4. 如果不能联网, 可以打包导入再安装。
      5. 如果性能要求比较高,最好选用推荐的MPI实现和编译器。
      6. 自己编译的话,不同版本有差异,不一定通用的。
      posted in OpenFOAM
      程
      程迪
    • RE: OpenFOAM初级入门建议(2020年更新版)

      @hungryandfool 赞,语言层级越抽象,越容易搞清楚逻辑。也越容易把“概念一致性”搞明白。

      而且Python并不一定低效,实际的工作都调用线代库或者Kernel其实也挺快的。

      posted in OpenFOAM
      程
      程迪
    • RE: 虚拟机 效率

      docker 好像没有这么多问题。

      posted in OpenFOAM
      程
      程迪
    • RE: 最快速的上google的方式是什么?

      @夜阑烟寒
      NICE!

      posted in C斯达克
      程
      程迪
    • RE: 压力方程和压力修正方程

      @李东岳
      同问,到底谁是压力方程,Pressure Linked Equation(PLE), Pressure Poisson Equation(PPE)?

      PISO用的好像也不是PPE
      @litong189456
      另外OpenFOAM是Picard迭代,

      至少SIMPLE的PLE是玩的Discrete Then Split, DTS流程。投影法之类的玩的是STD流程。

      posted in OpenFOAM
      程
      程迪
    • RE: 压力方程和压力修正方程

      @litong189456
      UEqn是一个fvMatrix<T> 类,类成员里有个指针psi_,指向了U,所以UEqn.A, H计算都会引用U,U更新,A,H也更新。

      posted in OpenFOAM
      程
      程迪
    • RE: 压力方程和压力修正方程

      https://chengdi123000.github.io/2018/01/05/OpenFOAM的不可压缩流算法/#压力联系方程之到底是哪个方程

      posted in OpenFOAM
      程
      程迪