doubleの行列を​各セルに内包している​cell型変数の分解​方法

9 views (last 30 days)
Noruji Muto
Noruji Muto on 26 Feb 2021
Commented: Noruji Muto on 28 Feb 2021
お世話になります。
下記のコードにおいて、t_req1~3の各セルに内包されている1✕45の全行列をy軸、altitudeをx軸としたグラフを作ることを最終的な目標としているのですが、
t_req1~3の各セルをそれぞれ個別に取得する事ができずエラーが出て困っています。
どなたか良い解決方法をご存じないでしょうか?
needed_velocity_tというエクセルファイルを送付します。
計算の過程で、このエクセルファイルからの読み込みが必要です。
検証の際はご利用ください。
clear
close all
A=1;
vini=1530;
syms r
assume(r > 0);
row=3000;
m=(4/3)*pi*(r^3)*row;
ips=8.854*10^-12;
e=1.602*10^-19;
k=1.38*10^-23;
fai=[-100,-10,5,10,18];
Vg=5000;
d=0.08;
E=Vg/d ;
R=sqrt(A/pi);
D=2*R;
%i=1:1000:50000;
%0~50kmまでの必要高度
filename = 'needed_velocity_t.xlsx';%軌道の維持に必要な速度即ち第一宇宙速度
sheet = 1;
xlRange = 'C2:C52';
subsetA = xlsread(filename,sheet,xlRange);
nv1=subsetA;%必要速度[m/s]
%50~100kmまでの必要高度
filename = 'needed_velocity_t.xlsx';
sheet = 1;
xlRange = 'C52:C102';
subsetB = xlsread(filename,sheet,xlRange);
nv2=subsetB;%必要速度[m/s]
%200~250kmまでの必要高度
filename = 'needed_velocity_t.xlsx';
sheet = 1;
xlRange = 'C202:C252';
subsetC = xlsread(filename,sheet,xlRange);
nv3=subsetC;%必要速度[m/s]
%0-50kmまでの高度取得
filename = 'needed_velocity_t.xlsx';
sheet = 1;
xlRange = 'B2:B52';
subsetD = xlsread(filename,sheet,xlRange);
H1=subsetD;%高度[m]
%50-100kmまでの高度取得
filename = 'needed_velocity_t.xlsx';
sheet = 1;
xlRange = 'B52:B102';
subsetE = xlsread(filename,sheet,xlRange);
H2=subsetE;%高度[m]
%200-250kmまでの高度取得
filename = 'needed_velocity_t.xlsx';
sheet = 1;
xlRange = 'B202:B252';
subsetF = xlsread(filename,sheet,xlRange);
H3=subsetF;%高度[m]
nd1=0.0039;
for i=1:51; %計算中では単位はm/s
around1(i)=(((i-1)*1000)+(10921*10^3))*pi; %任意高度における一周の距離
end
t_around1= around1/vini;%初速のままの時一周にかかる時間
tin1=t_around1;
Vin=A*vini*tin1; %取り込み口が通る体積
nin1=Vin.'*nd1; %取り込み時間中に取り込む粒子数
C=4*pi*ips*r; %キャパシタンス\
q=fai*C; %粒子の電荷量
%ここから排気速度の話
FE=E*q; %電場から受ける力
%必要な周回数の算出
int_q=int(q,[0.3*10^-3 4.96*10^-3],'IgnoreAnalyticConstraints', true)
d_q=double(int_q);
dq=d_q;
dltv=sqrt((2*abs(dq)*E*d)/m)-vini;
% bunbo=nv1.^2
% m_n=(bunshi.')/bunbo
% m_n_trans=m_n(:,1)
mdot1_n=nv1*A*nd1;
%int_mdot1=int(mdot1,[0.3*10^-3 4.96*10^-3],'IgnoreAnalyticConstraints', true)
%d_mdot1=double(int_mdot1)
nin_n=mdot1_n/m;
%nin_n=nin_n(:,1)
Vin_n=nin_n/nd1; %取り込み時間中に取り込む粒子数
tin1_n=Vin_n/(A*vini); %取り込み口が通る体積
t_around1_n=tin1_n;
around1_n=t_around1_n*vini;%
int_around1_n=int(around1_n,[0.3*10^-3 4.96*10^-3],'IgnoreAnalyticConstraints', true)
d_around1_n=double(int_around1_n);
% gom=pinv(around1_n)
nc1=(d_around1_n.')./around1;
%nd2の部
nd2=0.0036;
for i=1:51;
around2(i)=(((i+49)*1000)+(10921*10^3))*pi; %任意高度における一周の距離
end
t_around2= around2/vini;
tin2=t_around2;
Vin=A*vini*tin2; %取り込み口が通る体積
nin2=Vin.'*nd2; %取り込み時間中に取り込む粒子数
% C=4*pi*ips*r %キャパシタンス\
% q=fai*C %粒子の電荷量
% %ここから排気速度の話
% FE=E*q %電場から受ける力
%必要な周回数の算出nd2
dltv=sqrt((2*abs(dq)*E*d)/m)-vini;
% bunbo=nv2.^2
% m_n=(bunshi.')/bunbo
% m_n_trans=m_n(:,1)
mdot2_n=nv2*A*nd2;
%int_mdot1=int(mdot1,[0.3*10^-3 4.96*10^-3],'IgnoreAnalyticConstraints', true)
%d_mdot1=double(int_mdot1)
nin2_n=mdot2_n/m;
%nin_n=nin_n(:,1)
Vin2_n=nin2_n/nd2; %取り込み時間中に取り込む粒子数
tin2_n=Vin2_n/(A*vini); %取り込み口が通る体積
t_around2_n=tin2_n;
around2_n=t_around2_n*vini;%
int_around2_n=int(around2_n,[0.3*10^-3 4.96*10^-3],'IgnoreAnalyticConstraints', true)
d_around2_n=double(int_around2_n);
% gom=pinv(around1_n)
nc2=(d_around2_n.')./around2;
%nd3の部
nd3=0.0018;
for i=1:51;
around3(i)=(((i+199)*1000)+(10921*10^3))*pi; %任意高度における一周の距離
end
t_around3= around3/vini;
tin3=t_around3;
Vin3_n=A*vini*tin3; %取り込み口が通る体積
nin3_n=Vin.'*nd3; %取り込み時間中に取り込む粒子数
% C=4*pi*ips*r %キャパシタンス\
% q=fai*C %粒子の電荷量
% %ここから排気速度の話
% FE=E*q %電場から受ける力
%必要な周回数の算出nd3
dltv=sqrt((2*abs(dq)*E*d)/m)-vini;
% bunbo=nv2.^2
% m_n=(bunshi.')/bunbo
% m_n_trans=m_n(:,1)
mdot3_n=nv3*A*nd3;
%int_mdot1=int(mdot1,[0.3*10^-3 4.96*10^-3],'IgnoreAnalyticConstraints', true)
%d_mdot1=double(int_mdot1)
nin3_n=mdot3_n/m;
%nin_n=nin_n(:,1)
Vin3_n=nin3_n/nd3; %取り込み時間中に取り込む粒子数
tin3_n=Vin3_n/(A*vini); %取り込み口が通る体積
t_around3_n=tin3_n;
around3_n=t_around3_n*vini;%
int_around3_n=int(around3_n,[0.3*10^-3 4.96*10^-3],'IgnoreAnalyticConstraints', true)
d_around3_n=double(int_around3_n);
% gom=pinv(around1_n)
nc3=(d_around3_n.')./around3;
%必要周回数を得るのに必要な周回時間
Md=[1,2,3,10,12];%乾燥重量
Md1=Md(1,1);
Md2=Md(1,2);
Md3=Md(1,3);
Md4=Md(1,4);
Md5=Md(1,5);
%ΔVを満足するのに必要な粒子重量
DV=[20,40,80,120];%ΔV
DV1=DV(1,1);
DV2=DV(1,2);
DV3=DV(1,3);
DV4=DV(1,4);
ar=[0.3:0.6:1:5:10:45]*(10^-6);%新たに仮定した粒子半径
Vex=[125.04 74.882 51.356 15.736 9.4846 3.1883];
Vex1=Vex(1,1);
Vex2=Vex(1,2);
Vex3=Vex(1,3);
Vex4=Vex(1,4);
Vex5=Vex(1,5);
Vex6=Vex(1,6);
count=1;
for i=1:1:4;
for j=1:6;
for k=1:5;
W=Md(k)*(exp(DV(i)/Vex(j))-1);%ΔVを満足するのに必要な粒子重量
n=W./((4/3)*(ar.^3)*row);%ΔVを満足するのに必要な粒子の個数
%%
%各Mdと各粒子数密度の組み合わせにおける必要周回数
n_req1=n./nin1;%各Mdと各粒子数密度の組み合わせにおける必要周回数
n_req2=n./nin2;
n_req3=n./nin3_n;
%%
%必要周回数を満たすのに必要な航行時間
t_req1=t_around1*n_req1;%0~50km
t_req2=t_around2*n_req2;%50~100km
t_req3=t_around3*n_req3;%200~250km
NUM(count,1:11)={i,j,Vex(j),W,n,n_req1,n_req2,n_req3,t_req1,t_req2,t_req3};
count=count+1;
end
end
end
count2=count-1;
%ΔVを満足するのに必要な粒子の個数
n=NUM(1:count2,5);
%各Mdと各粒子数密度の組み合わせにおける必要周回数
n_req1=NUM(1:count2,6);
n_req2=NUM(1:count2,7);
n_req3=NUM(1:count2,8);
%必要周回数を満たすのに必要な航行時間
t_req1=NUM(1:count2,9);%0~50km
t_req2=NUM(1:count2,10);%50~100km
t_req3=NUM(1:count2,11);%200~250km
%%
%高度を行列にしたもの
count3=1;
for i=0:1:49;
altitude1=1+i;%0~50,50~100km
alt1(count3,1)={altitude1};
count3=count3+1;
end
count3.5=1;
for i=49:1:99;
altitude1=1+i;50~100km
alt1(count3,1)={altitude1};
count3.5=count3.5+1;
end
count4=count3+1;
for i=199:1:249;%200~250km
altitude2=1+i;
alt2(count4,1)={altitude2};
count4=count4+1;
end
altitude=cat(1,alt1,alt2);
count5=1;
for i=1:1:45;
a(i)=t_req1(1,i);
b(1,i)=cell2mat(a(i));%%%%%%%エラー発生個所%%%%%%%%%
%B(count5,45)=(b);
count5=count5+1;
end
現状で発生しているエラーは下記のとおりです
左辺のインデックスが右辺とサイズが適合しないため、代入は実行できません。
エラー: needed_cicle_vs_needed_velocity_ver6 (line 248)
b(1,i)=cell2mat(a(i));

Accepted Answer

Hernia Baby
Hernia Baby on 26 Feb 2021
Edited: Hernia Baby on 26 Feb 2021
エラーの原因はサイズの違いです。
i = 1で見てみると
size(cell2mat(a(i))) は 1×45に対し、size(b(1,i))は1×1になります。
以下のように差し替えてみてください。
count5=1;
for i=1:1:45
a(i)=t_req1(i,1); %ここも変更の必要あり
b(i,:)=cell2mat(a(i)); %%%%%%%エラー発生個所%%%%%%%%%
%B(count5,45)=(b);
count5=count5+1;
end
うまく回らなかったので、234~239行目も変えています。
count3p5=1; % 3.5だとエラー
for i=49:1:99
altitude1=1+i; % 50~100km コメントアウトし忘れ
alt1(count3,1)={altitude1};
count3p5=count3p5+1;
end
余談ですが、for文の条件はセミコロンでしめなくても大丈夫です。
  2 Comments
Hernia Baby
Hernia Baby on 26 Feb 2021
補足
252行目の
a(i)=t_req1(i,1);
にする理由
【理由】
t_req1が120×1の行列なので
i = 2のとき、t_req1(1,2)は存在しません
【他の方法】
他の方法として転置してあげるとうまくいきます
転置は以下の通り
t_req1 = t_req1';
Noruji Muto
Noruji Muto on 28 Feb 2021
ありがとうございます。
ためしに走らせてみたところ問題は解決されていることを確認しました!

Sign in to comment.

More Answers (0)

Community Treasure Hunt

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

Start Hunting!