2

我想对一个 3D 函数进行矢量化,但该函数没有解析表达式。例如,我可以对函数进行矢量化

 F(x, y, z) = (sin(y)*x, z*y, x*y)

通过做类似的事情

 function out = Vect_fn(x, y,z)
     out(1) = x.*sin(y);
     out(2) = z.*y;
     out(3) = x.*y;
 end

然后运行脚本

 a = linspace(0,1,10);
 [xx, yy, zz] = meshgrid(a, a, a);
 D = Vect_fn(xx, yy, zz)

但是,假设函数没有解析表达式,例如

 function y = Vect_Nexplicit(y0)
      %%%%%%y0 is a 3x1 vector%%%%%%%%%%%%%%
      t0 = 0.0; 
      tf = 3.0;
      [t, z] = ode45('ODE_fn', [t0,tf], y0);
      sz = size(z);
      n = sz(1);
      y = z(n, :);
 end

whereODE_fn只是一些在 ODE 右侧吐出的函数。因此,该函数仅求解 ODE,因此该函数不明确。当然,我可以使用for循环,但速度较慢(尤其是在 Octave 中,我更喜欢它,因为它lsode可以解决 ODE)

尝试类似的东西

a = linspace(0,1,10);
 [xx, yy, zz] = meshgrid(a, a, a);
 D = Vect_Nexplicit(xx, yy, zz)

不起作用。这里也是 ODF_fn 的代码:

 function ydot = ODE_fn(t, yin)

 A = sqrt(3.0);
 B = sqrt(2.0);
 C = 1.0;

 x = yin(1, 1);
 y = yin(2,1);
 z = yin(3, 1);

 M = reshape(yin(4:12), 3, 3);

 ydot(1,1) = A*sin(yin(3)) + C*cos(yin(2));
 ydot(2,1) = B*sin(yin(1)) + A*cos(yin(3));
 ydot(3,1) = C*sin(yin(2)) + B*cos(yin(1));

 DV = [0 -C*sin(y) A*cos(z); B*cos(x) 0 -A*sin(z); -B*sin(x)      C*cos(y) 0];

 Mdot = DV*M;

 ydot(4:12,1) = reshape(Mdot, 9, 1);


 end
4

1 回答 1

0

您可以使用 ode45 求解微分方程组,因此如果 ODE_fn 被矢量化,您可以使用您的方法。y0 也只需要是一个向量。

您可以创建一个 y0,即 [x1, ..., xn, y1, ..., yn, z1, ..., zn, M1_1-9, ... ,Mn_1-9] 然后用于 x, y, z 只是适当的索引,即 1:n、n+1:2*n、2*n+1:3n。然后使用 reshape(yin(3*n+1:end),3,3,n)。但我不确定如何向量化矩阵乘法。

于 2012-08-09T10:49:53.797 回答