如何在openfoam中求解一元高次方程

星星星星晴

如图
k1 k2 $We_0$ 都是常数,我想求解 $r_bu$

原文中是这么说的
The solution was numerically calculated using the linear interpolation method.

我使用了二分法求该方程，但是问题是现在不知道该如何将其放入CloudFunction中。
不知道哪里去define 一些东西。

星星星星晴

#include <iostream>
#include "fvCFD.H"
#include <stdio.h>
#include <math.h>

using namespace std;

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

long double f(float x){
	long double f;
	f= 0.75/1.414*8.5*0.45*5.0 * pow(x,3) + pow(x,2) - 1.0;
	return f;
}




// Main program:

int main()
{

long double x1,x2,x0,t,counter = 0.0;
x1 = -10;
x2 = 10;
do{
	x0=(x1+x2)/2.0;
	if( f(x0) * f(x1) > 0)
	{
		x1=x0;
	}
	else if (f(x0) * f(x1) < 0)
	{
		x2 = x0;
	}
t =fabs( f(x0));
counter +=1;
cout <<"Counter = "<< counter <<nl;
cout <<"t = "<< t <<nl;
if (counter > 100000)
{
	break;
}
}while(t> 1e-5);

cout<<"x0 = "<<x0<<nl;

cout << " solution = " << x0<<nl;
cout << "Hello World "<<nl;
return 0;

}


// ************************************************************************* //

星星星星晴

或者说有没有其他的解法？
例如该方法能不能解？
该如何修改该部分code呢？？

https://www.cfd-online.com/Forums/openfoam-programming-development/173789-how-solve-nonlinear-equation.html

麻烦大佬指点谢谢！

东岳

你的方程就是求$x^{3.5}+x^2-1=0$的根吧？

但是问题是现在不知道该如何将其放入CloudFunction中。

不知道哪里去define 一些东西。

你是要放在某个C程序里么还是什么地方？

星星星星晴

放在lagrangian 下面的submodel中

东岳

求根可以直接写成任意一个rootFinding函数，可以放在任何subModels下面的.C文件里面？

星星星星晴

@东岳大哥，这对我来说有点超纲

星星星星晴

求解该方程会在下面的这个class里面，但是我不知道在哪能定义f那个函数

template<class CloudType>
void Foam::IVT41_StochasticDDCollModel<CloudType>::postMove
(
    parcelType& p,
    const label cellI,
    const scalar dt,
    const point&,
    bool&
) {

//*code

}

星星星星晴

@东岳
上图的函数会在lagrangian/intermediate/submodel/CloudFunctionObjects 里，是自定义的一个模型

星星星星晴

@东岳是我着相了。。。我好傻，直接算就好了，干嘛搞一个class function。。。

scalar k1 = 8.5;
scalar k2 = 0.45;
scalar rbu_avg;
scalar x0 = 0.0;
scalar x2_rlarge = Dmax/2.0;
scalar x1_rsmall = Dmin/2.0;
scalar sat_counter = 0;
scalar eps = 0.0;
do
{
	scalar x0 = (x2_rlarge + x1_rsmall) / 2.0;
	scalar x0_ytempavg   = 0.75/sqrt(x0) * k1 *k2 * We0 * pow(x0,7.0/2.0) + pow(x0,2.0) - 1.0 ;
	scalar x2_ytemplarge = 0.75/sqrt(x2_rlarge) * k1 *k2 * We0 * pow(x2_rlarge,7.0/2.0) + pow(x2_rlarge,2.0) - 1.0 ;
	scalar x1_ytempsmall = 0.75/sqrt(x1_rsmall) * k1 *k2 * We0 * pow(x1_rsmall,7.0/2.0) + pow(x1_rsmall,2.0) - 1.0 ;

	if( x0_ytempavg * x1_ytempsmall > 0)
		{
			x1_rsmall=x0;
		}
		else if (x0_ytempavg * x1_ytempsmall < 0)
		{
			x2_rlarge = x0;
		}
	eps = fabs( x0_ytempavg);

	sat_counter +=1;

	if (sat_counter > 100000)
	{
		break;
	}
}while(eps > 1e-6);

rbu_avg = x0;

东岳

@星星星星晴是的