Problem using quad2d with vector

1 visualizzazione (ultimi 30 giorni)
Scott Lines
Scott Lines il 26 Ott 2016
Risposto: Torsten il 26 Ott 2016
Using below code with the variable t doesn't work, matrix dimensions don't agree, if the variable is replaced with a single number it does work (though that's not particularly useful), I've tried replacing it with a loop to no avail.
What am i missing?
lm = 1.75;
Cm = 3.2E6;
fe = 1380;
H = 13.5;
% time and space
y = 0.0;
ro = 0.1;
x = ro;
ts = [0.1 97/24];
Mt = 100;
% t = logspace(log10(ts(1)*86400),log10(ts(2)*86400),Mt); % [s]
t = linspace(ts(1)*86400, ts(2)*86400, Mt); % [s]
% calculated parameters
QL = fe/H;
z = H/2;
Dt = lm/Cm;
r = sqrt(x.^2+y.^2);
T = QL*(1/4/pi/lm)*...
(quad2d(@(f,ze) erfc(sqrt(r^2+ro^2-2*r*ro*cos(f)+(z-ze).^2)...
/2/sqrt(Dt*t))/pi./sqrt(r^2+ro^2-2*r*ro*cos(f)+(z-ze).^2),0,pi,0,H)-...
quad2d(@(f,ze) erfc(sqrt(r^2+ro^2-2*r*ro*cos(f)+(z-ze).^2)...
/2/sqrt(Dt*t))/pi./sqrt(r^2+ro^2-2*r*ro*cos(f)+(z-ze).^2),0,pi,-H,0));

Risposta accettata

KSSV
KSSV il 26 Ott 2016
You have to take care of element wise operations in that big expression T. As T is very big and clumsy, I have modified it to a loop, as below.
clc; clear all ;
lm = 1.75;
Cm = 3.2E6;
fe = 1380;
H = 13.5;
% time and space
y = 0.0;
ro = 0.1;
x = ro;
ts = [0.1 97/24];
Mt = 100;
% t = logspace(log10(ts(1)*86400),log10(ts(2)*86400),Mt); % [s]
tt = linspace(ts(1)*86400, ts(2)*86400, Mt); % [s]
% calculated parameters
QL = fe/H;
z = H/2;
Dt = lm/Cm;
r = sqrt(x.^2+y.^2);
T = zeros(size(tt)) ;
for i = 1:length(tt)
t = tt(i) ;
T(i) = QL*(1/4/pi/lm)*...
(quad2d(@(f,ze) erfc(sqrt(r^2+ro^2-2*r*ro*cos(f)+(z-ze).^2)...
/2/sqrt(Dt*t))/pi./sqrt(r^2+ro^2-2*r*ro*cos(f)+(z-ze).^2),0,pi,0,H)-...
quad2d(@(f,ze) erfc(sqrt(r^2+ro^2-2*r*ro*cos(f)+(z-ze).^2)...
/2/sqrt(Dt*t))/pi./sqrt(r^2+ro^2-2*r*ro*cos(f)+(z-ze).^2),0,pi,-H,0));
end
Though it is working, it is throwing some warning. I suggest you to check the expression for T.
  1 Commento
Scott Lines
Scott Lines il 26 Ott 2016
Thank you, I think I will need to go over the expression again but that's solved the first problem.

Accedi per commentare.

Più risposte (1)

Torsten
Torsten il 26 Ott 2016
There is no other way than to loop over the elements of t and call quad2d for each element separately.
MATLAB's "integral" has the option "ArrayValued=True", but as far as I know, integral2 or quad2d haven't.
Best wishes
Torsten.

Tag

Community Treasure Hunt

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

Start Hunting!

Translated by