How to convert into double?

2 次查看（过去 30 天）
Max
Max 2016-2-5
评论： Walter Roberson 2016-2-6
Hello,
I wrote the following code:
tfail = [5571.760,5573.742,5654.457,6079.693,6081.927,6172.915,6515.064,6517.515,6617.308,7095.558,7098.298,7209.831,7530.929,7533.885,7654.224,7966.300,7969.472,8098.617,8401.671,8405.059,8543.009,8982.166,8985.843,9135.533,9852.908,9857.017,10024.38,10868.774,10873.387,11061.234];
n=length(tfail);
beta_hat = 4.2915822;
B_hat = 1861.6186657;
C_hat = 58.9848692;
syms t B beta C
y(t) = (exp(-B/((heaviside(t)-heaviside(t-2000))*(330)+(heaviside(t-2000)-heaviside(t-3000))*(350)+...
(heaviside(t-3000)-heaviside(t-14000))*(390))))/C;
ogL=0;
for i=1:n
tfail(i);
I(i) = int(y(t),t,0,tfail(i));
y_new(i)=subs(y,t,tfail(i));
logL =logL+log((beta*y_new(i)*(I(i))^(beta-1))*exp(-((I(i))^beta)));
end
p = int(y(t),t,0,14000);
u = beta*log(p);
du_dB = diff(u,B);
du_dbeta = diff(u,beta);
du_dC = diff(u,C);
du_dB_sub = subs(du_dB,{beta,B,C},{beta_hat,B_hat,C_hat});
du_dbeta_sub = subs(du_dbeta,{B,C},{B_hat,C_hat});
du_dC_sub = subs(du_dC,{beta,B,C},{beta_hat,B_hat,C_hat});
v=[beta;B;C];
H=hessian(logL,v);
H_negatv=-1*H;
now I would like to calculate the inverse of H_negatv by using:
H_inverse=inv(H_negatv);
But that doesn´t work. So I tried out:
h = 1\H_negatv.
That´s good so far.
But now I would do sth. like that:
w=subs(h,[beta,B,C],[beta_hat,B_hat,C_hat]);
F_direct = w;
In according to calculate:
Var_B_hat_direct = double(F_direct(2,2));
But I can´t do that in MATLAB.
Does somebody have an idea how to solve that problem?
  7 个评论
显示 5更早的评论
Max
Max 2016-2-6
编辑：Walter Roberson 2016-2-6
Dear John,
first of all, you´re right concerning the matlab command: 1\A. That doesn´t calculate the inverse matrix.
But instead of doing Ainv = inv(A), I also can do: A_inv = A^(-1):
A = [ 1 2 3; 2 3 4; 1 2 5]
A =
1 2 3
2 3 4
1 2 5
inv(A)
ans =
-3.5000 2.0000 0.5000
3.0000 -1.0000 -1.0000
-0.5000 0 0.5000
A_inv
ans =
-3.5000 2.0000 0.5000
3.0000 -1.0000 -1.0000
-0.5000 0 0.5000
My problem is that I can´t do it in MATLAB because of the "out of memory"!!!
Walter Roberson
Walter Roberson 2016-2-6
You are not looping substituting different values for the symbols in H_negatv, so what you should be doing is substituting in the numeric values and double() the result before doing the inverse, so that you are taking the inverse of a numeric 3 x 3 instead of a large symbolic matrix.

请先登录，再进行评论。

请先登录，再回答此问题。

采纳的回答

Walter Roberson
Walter Roberson 2016-2-5
Working notes for me.
rational = @(V) sym(V, 'r');
if ismember( exist('hessian'), [2, 3, 5, 6, 8]) %is it an executable function?
Hess = @(M,V) hessian(M,V)
else
Hess = @(M,V) maple('Student[VectorCalculus][Hessian]', M, maple('convert', V, 'list'));
end
tfail = rational([5571.760, 5573.742, 5654.457, 6079.693, 6081.927, 6172.915, 6515.064, 6517.515, 6617.308, 7095.558, 7098.298, 7209.831, 7530.929, 7533.885, 7654.224, 7966.300, 7969.472, 8098.617, 8401.671, 8405.059, 8543.009, 8982.166, 8985.843, 9135.533, 9852.908, 9857.017, 10024.38, 10868.774, 10873.387, 11061.234]);
n = length(tfail);
beta_hat = rational(4.2915822);
B_hat = rational(1861.6186657);
C_hat = rational(58.9848692);
syms t B beta C
y = (exp(-B/((heaviside(t) - heaviside(t-2000)) * (330) + (heaviside(t-2000) - heaviside(t-3000)) * (350) + (heaviside(t-3000) - heaviside(t-14000)) * (390))))/C;
Z = rational(0);
LogL = Z;
I = rational(zeros(1,n));
y_new = rational(zeros(1,n));
new_term = rational(zeros(1,n));
for i = 1:n
I(i) = simplify( int(y, t, Z, tfail(i)) );
y_new(i) = simplify( subs(y, t, tfail(i)) );
new_term(i) = log((beta * y_new(i) * (I(i))^(beta-1)) * exp(-((I(i))^beta)));
logL = logL + new_term(i);
end
p = int(y, t, Z, rational(14000));
u = beta * log(p);
du_dB = diff(u, B);
du_dbeta = diff(u, beta);
du_dC = diff(u, C);
du_dB_sub = subs(du_dB, {beta, B, C}, {beta_hat, B_hat, C_hat});
du_dbeta_sub = subs(du_dbeta, {B,C}, {B_hat,C_hat});
du_dC_sub = subs(du_dC, {beta,B,C}, {beta_hat,B_hat,C_hat});
v = [beta; B; C];
H = Hess(logL, v);
%the next will probably fail, running out of memory
Hs = simplify(H);
H_negatv = -1*Hs;
  1 个评论
Walter Roberson
Walter Roberson 2016-2-6
The above is not intended as a solution: it is a recoding for compatibility with the Maple interface. On my system it will run out of memory attempting to simplify H. And if you skip the simplify() step and go on to taking the inv() it will run out of memory trying to take the inverse.
Your expression is simply too large to take a reasonable symbolic inverse of. But since you have specific numeric values to substitute in, you should substitute those specific numeric values into H, take double() of the result (since it will be completely numeric) and inv() that purely numeric interface.
Except, of course, we do not recommend that you take inv() for any reason other than proving that you can.

请先登录，再进行评论。

更多回答（0 个）

类别

在 Help Center 和 File Exchange 中查找有关 Linear Algebra 的更多信息

标签

产品

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by