力学分析开源软件简介.docx

Palabos (/)PalaBos的是一款高效的流体模拟及其建模库，开发基于C的STL(标准模板库)，有极强的拓展性！尽管其源代码是开放，但是基于PalaBos的FlowKit公司已于2011年9月开始运营（/），主要为流体力学相关领域提供解决方案，并定制软件。主要的开发者为我的日内瓦朋友Jonas Latt博士，另外一个重要开发成员Orestis博士也是我的合作者和好朋友，其主要的贡献在于湍流模型和多块加密的代码的开发。在版本1.0中，目前二维的多块加密是可用的，三维的曲面边界可用，需要提供stl几何文件（参：examples/showCases/aneurysm）。PalaBos的主要特点在于，其在并行结构上采取并行机制与模型分离的方式，使得应用建模与并行机制不相关。这也使得PalaBos的易于扩展。下面举例来说明其代码特点：对于二维计算下面两个基本的文件必须包括#include palabos2D.h#include palabos2D.hh#include #include #include #include #include 基本的名字空间using namespace plb;using namespace plb:descriptors;using namespace std;typedef double T;基本模型的描述。对于PalaBos，众多模型的应用，都是通过DnQmDescriptor来描述的。用户可自定义！ #define DESCRIPTOR D2Q9Descriptor初场的定义，建议使用这种方法T poiseuilleVelocity(plint iY, IncomprFlowParam const& parameters) T y = (T)iY / parameters.getResolution();return 4.*parameters.getLatticeU() * (y-y*y);压力的定义T poiseuillePressure(plint iX, IncomprFlowParam const& parameters) T Lx = parameters.getNx()-1;T Ly = parameters.getNy()-1;return 8.*parameters.getLatticeNu()*parameters.getLatticeU() / (Ly*Ly) * (Lx/(T)2-(T)iX);密度场的定义T poiseuilleDensity(plint iX, IncomprFlowParam const& parameters) return poiseuillePressure(iX,parameters)*DESCRIPTOR:invCs2 + (T)1;根据坐标进行速度初始化templateclass PoiseuilleVelocity public:PoiseuilleVelocity(IncomprFlowParam parameters_): parameters(parameters_) void operator()(plint iX, plint iY, Array& u) const u0 = poiseuilleVelocity(iY, parameters);u1 = T();private:IncomprFlowParam parameters;根据坐标进行密度初始化templateclass PoiseuilleDensity public:PoiseuilleDensity(IncomprFlowParam parameters_): parameters(parameters_) T operator()(plint iX, plint iY) const return poiseuilleDensity(iX,parameters);private:IncomprFlowParam parameters;零速度场初始化templateclass PoiseuilleDensityAndZeroVelocity public:PoiseuilleDensityAndZeroVelocity(IncomprFlowParam parameters_): parameters(parameters_) void operator()(plint iX, plint iY, T& rho, Array& u) const rho = poiseuilleDensity(iX,parameters);u0 = T();u1 = T();private:IncomprFlowParam parameters;枚举类型，表示边界类型enum InletOutletT pressure, velocity;建立几何void channelSetup( MultiBlockLattice2D& lattice,IncomprFlowParam const& parameters,OnLatticeBoundaryCondition2D& boundaryCondition,InletOutletT inletOutlet )const plint nx = parameters.getNx();const plint ny = parameters.getNy();下边界boundaryCondition.setVelocityConditionOnBlockBoundaries (lattice, Box2D(0, nx-1, 0, 0) );上边界boundaryCondition.setVelocityConditionOnBlockBoundaries (lattice, Box2D(0, nx-1, ny-1, ny-1) );设置相关边界if (inletOutlet = pressure) boundaryCondition.setPressureConditionOnBlockBoundaries (lattice, Box2D(0,0, 1,ny-2) );boundaryCondition.setPressureConditionOnBlockBoundaries (lattice, Box2D(nx-1,nx-1, 1,ny-2) );else boundaryCondition.setVelocityConditionOnBlockBoundaries (lattice, Box2D(0,0, 1,ny-2) );boundaryCondition.setVelocityConditionOnBlockBoundaries (lattice, Box2D(nx-1,nx-1, 1,ny-2) );设置常密度和速度，边界的速度和密度同时被施加setBoundaryDensity (lattice, lattice.getBoundingBox(),PoiseuilleDensity(parameters) );setBoundaryVelocity (lattice, lattice.getBoundingBox(),PoiseuilleVelocity(parameters) );初始化平衡态initializeAtEquilibrium (lattice, lattice.getBoundingBox(),PoiseuilleDensityAndZeroVelocity(parameters) );初始化模拟系统lattice.initialize();图片文件保存void writeGif(MultiBlockLattice2D& lattice, plint iter)const plint imSize = 600;ImageWriter imageWriter(leeloo);imageWriter.writeScaledGif(createFileName(u, iter, 6),*computeVelocityNorm(lattice),imSize, imSize );VTK 保存void writeVTK(MultiBlockLattice2D& lattice,IncomprFlowParam const& parameters, plint iter)T dx = parameters.getDeltaX();T dt = parameters.getDeltaT();VtkImageOutput2D vtkOut(createFileName(vtk, iter, 6), dx);vtkOut.writeData(*computeVelocityNorm(lattice), velocityNorm, dx/dt);vtkOut.writeData(*computeVelocity(lattice), velocity, dx/dt);计算均方根误差T computeRMSerror ( MultiBlockLattice2D& lattice,IncomprFlowParam const& parameters )MultiTensorField2D analyticalVelocity(lattice);setToFunction( analyticalVelocity, analyticalVelocity.getBoundingBox(),PoiseuilleVelocity(parameters) );MultiTensorField2D numericalVelocity(lattice);computeVelocity(lattice, numericalVelocity, lattice.getBoundingBox();return 1./parameters.getLatticeU() *std:sqrt( computeAverage( *computeNormSqr(*subtract(analyticalVelocity, numericalVelocity) ) );主函数int main(int argc, char* argv) plbInit(&argc, &argv);global:directories().setOutputDir(./tmp/);基本参数，Orestis为我加入了另外的个构造函数，可直接在物理空间中定义这些参数！后续我将给出分析！IncomprFlowParam parameters(T) 2e-2, / uMax(T) 5., / Re60, / N3., / lx1. / ly);const T logT = (T)0.1;const T imSave = (T)0.5;const T vtkSave = (T)2.;const T maxT = (T)15.1;const InletOutletT inletOutlet = velocity;写日志文件writeLogFile(parameters, Poiseuille flow);使用MultiBlockLatticeMultiBlockLattice2D lattice (parameters.getNx(), parameters.getNy(),设置松弛参数new BGKdynamics(parameters.getOmega() );声明边界OnLatticeBoundaryCondition2D*boundaryCondition = createLocalBoundaryCondition2D();channelSetup(lattice, parameters, *boundaryCondition, inletOutlet);主循环for (plint iT=0; iT*parameters.getDeltaT()maxT; +iT) if (iT%parameters.nStep(imSave)=0) pcout Saving Gif . 0) pcout Saving VTK file . endl;writeVTK(lattice, parameters, iT);if (iT%parameters.nStep(logT)=0) pcout step iT ; t= iT*parameters.getDeltaT() ; RMS error= computeRMSerror(lattice, parameters);Array uCenter;lattice.get(parameters.getNx()/2,parameters.getNy()/2).computeVelocity(uCenter);pcout ; center velocity= uCenter0/parameters.getLatticeU() endl;碰撞迁移，此函数不再有bool类型参数，其不同于Openlblattice.collideAndStream();删除边界对象delete boundaryCondition;由此上面代码可以看出，PalaBos的结构是非常清晰的，用户在实现过程也异常的简单！以上资料出自/s/blog_6a40fba10100znv9.htmlParaview (/) Par