reactingFoam/reactingDNS计算一维球形火焰传播

[bany]

（One dimensional spherical flame initiation）
我们最近在用reactingFoam、reactingDNS（考虑了组分的详细输运特性）来计算1D球形火焰传播，并且想与ASURF对照。在该算例中，火焰是在固定区域内沉积能量来引发的，该能量沉积区域位于1D计算域的左侧，长度rig=0.2mm，持续时长tig=0.2ms。该能量直接加入能量方程：power=energy/(4/3X3.14Xrig^3Xtig)Xexp(-(r/rig)^6)（X表示*）。在算例中使用topoSet工具降能量沉积区域disZone圈出来，求解器内部读取该disZone中的网格变化，然后对能量源项进行赋值。

首先创建能量源项sparkE，读取算例中给的能量沉积半径dis_radius，持续时间dis_duration，能量dis_energy。
volScalarField sparkE
(
    IOobject
    (
        "sparkE",
        runTime.timeName(),
        mesh,
        IOobject::NO_READ,
        IOobject::AUTO_WRITE
    ),
    mesh,
	dimensionedScalar("sparkE", dimEnergy/dimTime/dimVolume, 0.0) // J/m3/s
);

IOdictionary chemistryDict
(
	IOobject
	(
		thermo.phasePropertyName("chemistryProperties"),
		thermo.db().time().constant(),
		thermo.db(),
		IOobject::MUST_READ,
		IOobject::NO_WRITE,
		false
	)
);

const dictionary& dischargeTypeDict = chemistryDict.subDict("sparkDischarge");

scalar dis_duration = dischargeTypeDict.lookupOrDefault<scalar>("dis_duration", 2E-4); // 200 us
scalar dis_radius = dischargeTypeDict.lookupOrDefault<scalar>("dis_radius", 2E-4); // 200 um
scalar dis_energy = dischargeTypeDict.lookupOrDefault<scalar>("dis_energy", 1E-3); // 1 mJ
//scalar dis_central = dischargeTypeDict.lookupOrDefault<scalar>("dis_central", 1E-3); // 1 mJ
scalar PI = 3.1415926;
scalar dis_paramTot = dis_energy/(4.0*PI/3.0*dis_radius*dis_radius*dis_radius*dis_duration);
//scalar dis_paramTot = dis_energy/(PI*Foam::sqrt(PI)*dis_radius*dis_radius*dis_radius*dis_duration);

在求解器主程序中，对能量赋值
	const label& disZoneID = mesh.cellZones().findZoneID("disZone");
	const labelList& cZone = mesh.cellZones()[disZoneID];
	forAll(mesh.C(),cellI)
	{
		sparkE[cellI] = 0.0;
	}
		
	if (disZoneID != -1 && cZone.size() > 0 && runTime.value() <= dis_duration)
	{
		forAll(cZone,cellI)
		{
			sparkE[cZone[cellI]] = dis_paramTot*std::exp(-PI/4.0*pow(mesh.C()[cZone[cellI]].x()/dis_radius,6));
//			sparkE[cZone[cellI]] = dis_paramTot*std::exp(-pow(mesh.C()[cZone[cellI]].x()/dis_radius,2));
		}
	}
	if ( (runTime.value()+runTime.deltaTValue()) > dis_duration && runTime.value() < dis_duration )
	{
		runTime.setDeltaT((dis_duration-runTime.value()));
	}
	Info << " sparkE max/min: "<<  max(sparkE).value() <<" "<< min(sparkE).value() << endl;

但是如果采用ASURF中不能点火的能量，在openfoam中竟然能点燃。
并且火焰半径（即火焰位置，定义为最大温度梯度）随时间的变化与ASURF完全不一样，图片左侧为openfoam计算结果，右侧为ASURF计算结果。请大家指点

算例可见附件：
1Dflame_case.zip