decomposePar运行时出现问题
-
不知道有没有哪位遇到过这样的情况,我在decomposePar的时候,出现这样的错误提示:
Processor 0: field transfer terminate called after throwing an instance of 'std::logic_error' _what(): basic_string::_S_consturct null not valid
入口边界条件是我自己定义的,功能是读取一个文件数据,进行运算之后赋值给这个边界条件的私有变量。我用的是simple分块,不破坏这个入口。哪里出错了呢可以给点建议吗?
PS:单核验证可以计算
-
可能我描述的问题不够准确,下面附上我的部分代码:
在.H文件中,定义了一些私有变量与私有变量的接口:// private perameters //- freestream Ma scalar Ma_; //- freestream T scalar Tinf_; scalar gamma_; scalar Prt_; //- munber of modes scalar nmodes_; //- smallest wavenumber scalar dxmin_; //- ratio of ke and kmin (in wavenumber) scalar wew1fct_; //- scalar amp_; //- kinetic viscosity of the smallest eddy scalar visc_; //- delta Inlet scalar deltaInlet_; //- freestream Velocity scalar Uinf_; //- freestream pressure scalar pInf_; //- read file IFstream IF_; //- rows of file data label ny_; //- MeanField vectorField Ubar_; //- primeField vectorField UPrime_; //- ReynoldStress tensor symmTensorField ReynoldStress_; label curTimeIndex_; // member function void calculateProfiles(IFstream&, vectorField&, symmTensorField&, scalar&); // Access scalar nmodes() const { return nmodes_; } scalar deltaInlet() const { return deltaInlet_; } scalar& deltaInlet() { return deltaInlet_; } scalar amp() const { return amp_; } scalar dxmin() const { return dxmin_; } scalar visc() const { return visc_; } scalar wew1fct() const { return wew1fct_; } scalar Ma() const { return Ma_; } scalar Tinf() const { return Tinf_; } scalar gamma() const { return gamma_; } scalar Prt() const { return Prt_; } scalar Uinf() const { return Uinf_; } scalar pInf() const { return pInf_; } IFstream IF() const { return IF_; } IFstream& IF() { return IF_; } vectorField& Ubar() { return Ubar_; } vectorField Ubar() const { return Ubar_; } symmTensorField& ReynoldStress() { return ReynoldStress_; } symmTensorField ReynoldStress() const { return ReynoldStress_; } label ny() { return ny_; }
然后在.C中,我的构造函数:
// constrctor Foam::syntheticVelocityFixedValueFvPatchField::syntheticVelocityFixedValueFvPatchField ( const fvPatch& p, const DimensionedField<vector, volMesh>& iF, const dictionary& dict ) : fixedValueFvPatchVectorField(p, iF), Ma_(readScalar(dict.lookup("Ma"))), Tinf_(readScalar(dict.lookup("Tinf"))), gamma_(dict.lookupOrDefault("gamma", 1.4)), Prt_(dict.lookupOrDefault("Prt", 0.9)), nmodes_(readScalar(dict.lookup("nmodes"))), dxmin_(readScalar(dict.lookup("dxmin"))), wew1fct_(dict.lookupOrDefault("wew1fct", 2.0)), amp_(dict.lookupOrDefault("amp", 1.452762113)), visc_(dict.lookupOrDefault("visc", 15.2e-6)), Uinf_(readScalar(dict.lookup("Uinf"))), pInf_(readScalar(dict.lookup("pInf"))), IF_(dict.lookup("pathName")), // stored meanStreamVelocity and ReynoldStress profile data ny_(dict.lookupOrDefault("ny", 535)), curTimeIndex_(-1) { if (dict.found("Ubar") && dict.found("ReynoldStress") && dict.found("deltaInlet")) { deltaInlet_= scalar(readScalar(dict.lookup("deltaInlet"))); Ubar_ = vectorField("Ubar", dict, p.size()); ReynoldStress_ = symmTensorField("ReynoldStress", dict, p.size()); } else { calculateProfiles(IF(), Ubar(), ReynoldStress(), deltaInlet()); } if (dict.found("UPrime")) { UPrime_ = vectorField("UPrime", dict, p.size()); } else { UPrime_ = vectorField(p.size(), vector::zero); } if (dict.found("value")) { fvPatchField<vector>::operator= ( vectorField("value", dict, p.size()) ); } else { fvPatchField<vector>::operator=(Ubar()); } }
这里调用的成员函数calculateProperties如下:
// member function void Foam::syntheticVelocityFixedValueFvPatchField::calculateProfiles ( IFstream& IF, vectorField& Ubar, symmTensorField& ReynoldStress, scalar& deltaInlet ) { if (!IF) { FatalErrorIn("void Foam::syntheticVelocityFixedValueFvPatchField::calculateMeanVelocity") << "can not find " << IF.name() << exit(FatalError); } else { scalarField UvdPlus(ny()), vPlus(ny()), wPlus(ny()), Ux(ny()), Uy(ny()), Uz(ny()); symmTensorField stressPlus(ny(), symmTensor::zero); scalarField yPlus(ny()), y(ny()); scalar yInf, ySup; scalar R(287.0), max_UvdPlus(26.457427), As(1.512e-06), Ts(120.0); scalar r=pow(Prt(), 1.0/3.0); //- calculate nuwInlet and Utao scalar Taw = Tinf()*(1.0+(gamma()-1.0)/2.0*r*sqr(Ma())); scalar Tw = (Taw - Tinf())/r + Tinf(); scalar rhow = pInf()/(R*Tw); scalar muw = As*sqrt(Tw)/(1.0+Ts/Tw); scalar nuwInlet = muw/rhow; scalar A = sqrt((gamma()-1.0)/2.0*Prt()*sqr(Ma())*Tinf()/Tw); scalar B = (1.0+sqrt(Prt())*(gamma()-1.0)/2.0*sqr(Ma()))*Tinf()/Tw - 1.0; scalar Utao = Uinf()/A*( asin((2.0*sqr(A)-B)/sqrt(sqr(B)+4.0*sqr(A))) + asin(B/sqrt(sqr(B)+4.0*sqr(A))) )/max_UvdPlus; // max(Uvd+) = 25.6803733; const scalarField T = this->patch().lookupPatchField<volScalarField, scalar>("T"); //- flag to decide whether to calculate deltaInlet bool delta=true; forAll(UvdPlus, patchI) { IF >> yPlus[patchI] >> UvdPlus[patchI] >> vPlus[patchI] >> wPlus[patchI] >> stressPlus[patchI].xx() >> stressPlus[patchI].yy() >> stressPlus[patchI].zz() >> stressPlus[patchI].xy(); // calculate dimensional y y[patchI] = yPlus[patchI]*nuwInlet/Utao; // calculate dimensional Ux scalar tmp = A*(UvdPlus[patchI]*Utao)/Uinf() - asin(B/sqrt(sqr(B)+4.0*sqr(A))); Ux[patchI] = (sin(tmp)*sqrt(sqr(B)+4.0*sqr(A)) + B)/(2.0*sqr(A))*Uinf(); Uy[patchI] = vPlus[patchI]*Utao; Uz[patchI] = wPlus[patchI]*Utao; if (Ux[patchI] > 0.99*Uinf() && delta) { deltaInlet = y[patchI]; delta = false; // deltaInlet calculated only once; } } vectorField Cf = this->patch().Cf(); Ubar = vectorField(Cf.size()); symmTensorField Reynold(Cf.size(), symmTensor::zero); ReynoldStress = symmTensorField(Cf.size(), symmTensor::zero); // interpolate at the cell forAll(Cf, faceI) { for(label patchI=0; patchI < y.size()-1; patchI++) { yInf = Cf[faceI][1]-y[patchI]; ySup = Cf[faceI][1]-y[patchI+1]; if(yInf*ySup <= 0) { Ubar[faceI][0] = yInf/(yInf-ySup)*Ux[patchI+1] - ySup/(yInf-ySup)*Ux[patchI]; Ubar[faceI][1] = yInf/(yInf-ySup)*Uy[patchI+1] - ySup/(yInf-ySup)*Uy[patchI]; Ubar[faceI][2] = yInf/(yInf-ySup)*Uz[patchI+1] - ySup/(yInf-ySup)*Uz[patchI]; Reynold[faceI] = yInf/(yInf-ySup)*stressPlus[patchI+1] - ySup/(yInf-ySup)*stressPlus[patchI]; // calculate dimensional Reynold stress scalar ksi = sqrt(Tw/T[faceI]); ReynoldStress[faceI].xx() = sqr(Reynold[faceI].xx()/ksi)*Utao; ReynoldStress[faceI].xy() = -sqr(Reynold[faceI].xy()/ksi)*Utao; ReynoldStress[faceI].yy() = sqr(Reynold[faceI].yy()/ksi)*Utao; ReynoldStress[faceI].zz() = sqr(Reynold[faceI].zz()/ksi)*Utao; break; } else // Cf.component(1) is higher than y.last() { Ubar[faceI][0] = Ux.last(); Ubar[faceI][1] = Uy.last(); Ubar[faceI][2] = Uz.last(); scalar ksi = sqrt(Tw/T[faceI]); ReynoldStress[faceI].xx() = sqr(stressPlus.last().xx()/ksi)*Utao; ReynoldStress[faceI].xy() = -sqr(stressPlus.last().xy()/ksi)*Utao; ReynoldStress[faceI].yy() = sqr(stressPlus.last().yy()/ksi)*Utao; ReynoldStress[faceI].zz() = sqr(stressPlus.last().zz()/ksi)*Utao; } } } scalar ymin = min(Cf.component(1)); scalar Uxmin = min(Ubar.component(0)); scalar dUxdy = Uxmin/ymin; scalar dUdy = sqr(Utao)/nuwInlet; Info << "dUxdy (from case simulation) =" << dUxdy << nl << "dUdy (from theory calculation) =" << dUdy << nl << "Utao =" << Utao << nl << "nuwInlet =" << nuwInlet << nl << "deltaInlet =" << deltaInlet << endl; } }
我在单核下测试时时可以计算的,但是想并行计算的时候,decomposePar过不去,终端输出信息如下:
/*---------------------------------------------------------------------------*\ | ========= | | | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox | | \\ / O peration | Version: 2.4.0 | | \\ / A nd | Web: www.OpenFOAM.org | | \\/ M anipulation | | \*---------------------------------------------------------------------------*/ Build : 2.4.0-dcea1e13ff76 Exec : decomposePar Date : Nov 13 2016 Time : 10:59:43 Host : "zwk" PID : 12485 Case : /home/Ubuntu/WORKPLACE4/SyntheticCase-test nProcs : 1 sigFpe : Enabling floating point exception trapping (FOAM_SIGFPE). fileModificationChecking : Monitoring run-time modified files using timeStampMaster allowSystemOperations : Allowing user-supplied system call operations // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // Create time Decomposing mesh region0 Create mesh Calculating distribution of cells Selecting decompositionMethod simple Finished decomposition in 0.04 s Calculating original mesh data Distributing cells to processors Distributing faces to processors Distributing points to processors Constructing processor meshes Processor 0 Number of cells = 23686 Number of faces shared with processor 1 = 1997 Number of faces shared with processor 1 = 58 Number of processor patches = 2 Number of processor faces = 2055 Number of boundary faces = 4243 Processor 1 Number of cells = 23686 Number of faces shared with processor 0 = 1997 Number of faces shared with processor 0 = 58 Number of faces shared with processor 2 = 1978 Number of faces shared with processor 2 = 52 Number of processor patches = 4 Number of processor faces = 4085 Number of boundary faces = 2189 Processor 2 Number of cells = 23686 Number of faces shared with processor 1 = 1978 Number of faces shared with processor 1 = 52 Number of faces shared with processor 3 = 1954 Number of faces shared with processor 3 = 11 Number of processor patches = 4 Number of processor faces = 3995 Number of boundary faces = 2329 Processor 3 Number of cells = 23686 Number of faces shared with processor 2 = 1954 Number of faces shared with processor 2 = 11 Number of faces shared with processor 4 = 1992 Number of faces shared with processor 4 = 57 Number of processor patches = 4 Number of processor faces = 4014 Number of boundary faces = 2276 Processor 4 Number of cells = 23686 Number of faces shared with processor 3 = 1992 Number of faces shared with processor 3 = 57 Number of faces shared with processor 5 = 1964 Number of faces shared with processor 5 = 43 Number of processor patches = 4 Number of processor faces = 4056 Number of boundary faces = 2208 Processor 5 Number of cells = 23686 Number of faces shared with processor 4 = 1964 Number of faces shared with processor 4 = 43 Number of faces shared with processor 6 = 1984 Number of faces shared with processor 6 = 41 Number of processor patches = 4 Number of processor faces = 4032 Number of boundary faces = 2330 Processor 6 Number of cells = 23685 Number of faces shared with processor 5 = 1984 Number of faces shared with processor 5 = 41 Number of faces shared with processor 7 = 1986 Number of faces shared with processor 7 = 55 Number of processor patches = 4 Number of processor faces = 4066 Number of boundary faces = 2218 Processor 7 Number of cells = 23685 Number of faces shared with processor 6 = 1986 Number of faces shared with processor 6 = 55 Number of processor patches = 2 Number of processor faces = 2041 Number of boundary faces = 4211 Number of processor faces = 14172 Max number of cells = 23686 (0.00105549% above average 23685.8) Max number of processor patches = 4 (14.2857% above average 3.5) Max number of faces between processors = 4085 (15.2978% above average 3543) Time = 0 dUxdy (from case simulation) =2.46405e+06 dUdy (from theory calculation) =2.46471e+06 Utao =23.1659 nuwInlet =0.000217736 deltaInlet =0.0144036 Processor 0: field transfer terminate called after throwing an instance of 'std::logic_error' what(): basic_string::_S_construct null not valid
还麻烦各位可以给个意见,这个问题就是这个边条造成的,但不知道如何解决,谢谢各位了
-
@Aeronastro 粗略看了一下,感觉问题出在
IF_
上。你的 dict 能贴上来吗? -
@wwzhao 好的,这个是我速度入口的边界条件
boundaryField { inlet { type syntheticVelocity; Ma 2.3; Tinf 104; Prt 0.9; nmodes 200; dxmin 0.01e-03; Uinf 552; pInf 3998.63; pathName "/home/Ubuntu/WORKPLACE4/syntheticCase-test/data"; value uniform (552 0 0); } other boundaries ...... }
由于刚开始计算的时候不知道Ubar,ReynoldStress和deltaInlet,在边条构造时需要通过那个calculateProfiles函数计算。
那个pathName 是我的文件所在的地方
其中的data文件中的数据是一列一列排布的,
0_1479171527178_data -
这个情况我倒是没有遇到过,不过我编了个synthetic的湍流生成条件,核比较少的时候patch像是一个完整的剖面,而decomposePar之后(没有遇到你这个问题)而且核数特别高(~500)之后,发现patch也被分割了的感觉。。我想说的是核数少的时候(~30)根本没发现patch上面会被不同的processor分割。
有一个计算圆形(是的,我的BC仅仅适用于圆形)中心的程序,我想不同processor算出来的ctr不同!
// Get range and orientation boundBox bb(patch().patch().localPoints(), true); vector ctr = 0.5*(bb.max() + bb.min());
看来得定义一个常量才稳妥
-
你好,我想请问,你怎么上传这些可以滑动的页面呀?