Assignment has more non-singleton rhs dimensions than non-singleton subscripts

1 visualización (últimos 30 días)
clear all clc
syms h;
v0=150:10:230;
c = [1 2 5 6 10 15 20 30 40];
th = zeros(length(v0), length(c));
tp = 0.1;
for i = 1:length(c)
for j = 1:length(v0)
t1 = (log(1+0.0509683996*c(i)*v0(j)))/c(i);
assume(h, 'real')
a = 1*10^10*c(i).^2;
b = 21962000000*c(i);
f = (-2*10^10)*c(i)*sqrt(3)-(4*10^8)*(v0(j).^2)-98100000*c(i)+36460944400;
d = (8*10^7)*(v0(j).^2)-(4*10^19)*sqrt(3)-157705560;
e = (4*10^6)*(v0(j).^2)+88400962361;
th (i,j) = vpasolve((a*h^4)+(b*h^3)+(f*h^2)+(d*h)-(e),h,[tp t1]);
t1 = 0;
end
end
  2 comentarios
Jan
Jan el 16 de Abr. de 2018
Editada: Jan el 16 de Abr. de 2018
Please post the complete error message, which mentions also, in which line the error occurs. This is better than letting the readers guess this important detail.
By the way, -2*10^10 is a multiplication and an expensive power operation, while -2e10 is a cheap constant. Compare the readability:
f = (-2*10^10)*c(i)*sqrt(3)-(4*10^8)*(v0(j).^2)-98100000*c(i)+36460944400;
f = -2e10 * c(i) * sqrt(3) - 4e8 * v0(j) .^ 2 - 981e5 * c(i) + 36460944400;
Muhammed AlSubhi
Muhammed AlSubhi el 16 de Abr. de 2018
Thank you for the hint, I just copied the code from maple and has no chance to edited the code.

Iniciar sesión para comentar.

Respuesta aceptada

Birdman
Birdman el 16 de Abr. de 2018

The output of vpasolve seems to be returning empty solutions. You can detect if the solutions are empty and assign an eps value instead but it seems all of your solutions return empty ones. That has to do something with your equation. Try this:

syms h;
v0=150:10:230;
c = [1 2 5 6 10 15 20 30 40];
th = zeros(length(v0), length(c));
tp = 0.1;
for i = 1:length(c)    
    for j = 1:length(v0)    
    t1 = (log(1+0.0509683996*c(i)*v0(j)))/c(i);    
    assume(h, 'real')    
    a = 1*10^10*c(i).^2;
    b = 21962000000*c(i);
    f = (-2*10^10)*c(i)*sqrt(3)-(4*10^8)*(v0(j).^2)-98100000*c(i)+36460944400;
    d = (8*10^7)*(v0(j).^2)-(4*10^19)*sqrt(3)-157705560;
    e = (4*10^6)*(v0(j).^2)+88400962361;
    res=vpasolve((a*h^4)+(b*h^3)+(f*h^2)+(d*h)-(e),h,[tp t1]);
    if isempty(res)
        th (i,j) = eps;%you need to assign a value, this is a very small one,eps=2.2e-16
    else
        th(i,j)=res;
    end
    t1 = 0; 
    end
end
  5 comentarios
Birdman
Birdman el 16 de Abr. de 2018
This piece of code differs mathematically from the initial code. Use the following code to get the results you want:
syms h;
v0=150:10:230;
c = [1 2 5 6 10 15 20 30 40];
th = zeros(length(v0), length(c));
tp = 0.1;
for i = 1:length(c)
for j = 1:length(v0)
t1 = (log(1+0.0509683996*c(i)*v0(j)))/c(i);
assume(h, 'real')
a = 1*10^10*c(i).^2;
b = 21962000000*c(i);
f = (-2*10^10)*c(i)*sqrt(3)-(4*10^8)*(v0(j).^2)-98100000*c(i)+36460944400;
d = (8*10^7)*(v0(j).^2)-(4*10^10)*sqrt(3)-157705560;
e = -(4*10^6)*(v0(j).^2)+88400962361;
res=vpasolve(a*h^4+b*h^3+f*h^2+d*h+e,h, [tp t1]);
if isempty(res)
th (i,j) = eps;%you need to assign a value, this is a very small one,eps=2.2e-16
else
th(i,j)=res;
end
t1 = 0;
end
end
Muhammed AlSubhi
Muhammed AlSubhi el 16 de Abr. de 2018
it did work thank you very much :) How can I plot c vo and th in a 3rd plot?

Iniciar sesión para comentar.

Más respuestas (0)

Etiquetas

Community Treasure Hunt

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

Start Hunting!

Translated by