用打靶法求解自然对流方程(热量和流量)

2024-06-01 11:59:39 发布

您现在位置:Python中文网/ 问答频道 /正文

TL;DR I've been implementing a python program to solve numerically equations for natural convection based on a particular similarity variable using runge-kutta 4 and the shooting method. However I don't get the right solutions when I plot it. Did I make a mistake somewhere ?

从自然对流的特殊情况出发,我们得到了这些相似方程。
第一个描述流体流动,第二个描述热流。
“Pr”是指Prandtl,它基本上是流体动力学中使用的无量纲数(Prandtl):

blasius problem for natural convection

这些方程受以下边界值的约束,因此板附近的温度大于边界层外的温度,并且流体速度远离边界层为0

boundary values

我一直在尝试用Runge Kutta 4和打靶法数值求解这些问题,将边值问题转化为初值问题。打靶法的实现方式是牛顿法

然而,我没有得到正确的解决方案。 正如你在下面所看到的,当我们离开极板时,温度(红色)正在升高,而它应该呈指数下降。 它对于流体速度(蓝色)更一致,但是我认为它应该上升得更快,然后下降得更快。这里的曲线更平滑

plot

事实上,我们有一个由2个耦合的ODE组成的系统。然而,现在,我只想找到两个初始值中的一个(例如,f''(0)=a,试图找到a),这样我们就有了一个边值问题的解决方案(shooting method)。一旦找到,我想我们就有了解决整个问题的办法

我想我应该管理这两个(f'(0)=a;θ(0)=b),但我不知道如何并行管理这两个。 最后想一提的是,如果我试图得到θ′(所以θ′(0))的初始值,我得不到正确的热量分布

代码如下:

"""
The goal is to resolve a 3rd order non-linear ODE for the blasius problem.
It's made of 2 equations (flow / heat)

f''' = 3ff'' - 2(f')^2 + theta
3 Pr f theta' + theta'' = 0

RK4 + Shooting Method
"""

import numpy as np
import math

from scipy.integrate import odeint
from scipy.optimize import newton

from edo_solver.plot import plot

from constants import PRECISION

def blasius_edo(y, t, prandtl):
  f = y[0:3]
  theta = y[3:5]
  return np.array([
    # flow edo
    f[1], # f' = df/dn
    f[2], # f'' = d^2f/dn^2
    - 3 * f[0] * f[2] + (2 * math.pow(f[1], 2)) - theta[0], # f''' = - 3ff'' + 2(f')^2 - theta,
    # heat edo
    theta[1], # theta' = dtheta/dn
    - 3 * prandtl * f[0] * theta[1], # theta'' = - 3 Pr f theta'
  ])

def rk4(eta_range, shoot):
  prandtl = 0.01

  # initial values
  f_init = [0, 0, shoot] # f(0), f'(0), f''(0)
  theta_init = [1, shoot] # theta(0), theta'(0)
  ci = f_init + theta_init # concatenate two ci

  # note: tuple with single argument must have "," at the end of the tuple
  return odeint(func=blasius_edo, y0=ci, t=eta_range, args=(prandtl,))

"""
if we have :
f'(t_0) = fprime_t0 ; f'(eta -> infty) = fprime_inf
we can transform it into :
f'(t_0) = fprime_t0 ; f''(t_0) = a

we define the function F(a) = f'(infty ; a) - fprime_inf
if F(a) has a root in "a",
then the solutions to the initial value problem with f''(t_0) = a
is also the solution the boundary problem with f'(eta -> infty) = fprime_inf

our goal is to find the root, we have the root...we have the solution.
it can be done with bissection method or newton method.
"""
def shooting(eta_range):
  # boundary value
  fprimeinf = 0 # f'(eta -> infty) = 0

  # initial guess
  # as far as I understand
  # it has to be the good guess
  # otherwise the result can be completely wrong
  initial_guess = 10 # guess for f''(0)

  # define our function to optimize
  # our goal is to take big eta because eta should approach infty
  # [-1, 1] : last row, second column => f'(eta_final) ~ f'(eta -> infty)
  fun = lambda initial_guess: rk4(eta_range, initial_guess)[-1, 1] - fprimeinf
  # newton method resolve the ODE system until eta_final
  # then adjust the shoot and resolve again until we have a correct shoot
  shoot = newton(func=fun, x0=initial_guess)

  # resolve our system of ODE with the good "a"
  y = rk4(eta_range, shoot)
  return y

def compute_blasius_edo(title, eta_final):
  ETA_0 = 0
  ETA_INTERVAL = 0.1
  ETA_FINAL = eta_final

  # default values
  title = title
  x_label = "$\eta$"
  y_label = "profil de vitesse $(f'(\eta))$ / profil de température $(\\theta)$"
  legends = ["$f'(\eta)$", "$\\theta$"]

  eta_range = np.arange(ETA_0, ETA_FINAL + ETA_INTERVAL, ETA_INTERVAL)

  # shoot
  y_set = shooting(eta_range)

  plot(eta_range, y_set, title, legends, x_label, y_label)

compute_blasius_edo(
  title="Convection naturelle - Solution de similitude",
  eta_final=10
)

Tags: thetoimportrangemethodinitialetawe
1条回答
网友
1楼 · 发布于 2024-06-01 11:59:39

您需要使用android.print包中的类

您需要使打印文档适配器类派生自PrintDocumentAdapter 这个类负责你打印的文档内容

要执行打印,请执行以下操作:

 // create and initalize your PrintDocumentAdapter
 ...

 // obtrain PrintManager instance
 PrintManager printManager =
     (PrintManager) context.getSystemService(Context.PRINT_SERVICE);

 // build PrintAttributes
 PrintAttributes.Builder builder = new PrintAttributes.Builder();
 builder. <set your print attributes such as: color/monochrome, media size, min margins>

 // print document 
 printManager.print(yourPrintJobName, yourPrintDocumentAdapter, builder.build());

注:安卓。打印要求API级别19及更高

相关问题 更多 >