Problem using quad2d with vector

1 vue (au cours des 30 derniers jours)
Scott Lines
Scott Lines le 26 Oct 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));

Réponse acceptée

KSSV
KSSV le 26 Oct 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 commentaire
Scott Lines
Scott Lines le 26 Oct 2016
Thank you, I think I will need to go over the expression again but that's solved the first problem.

Connectez-vous pour commenter.

Plus de réponses (1)

Torsten
Torsten le 26 Oct 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.

Catégories

En savoir plus sur Historical Contests dans Help Center et File Exchange

Tags

Community Treasure Hunt

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

Start Hunting!

Translated by