How can I plot the innermost contour?

7 views (last 30 days)
clc;
clear;
L = 100;
v_origin = 10;
v_opposite = 20;
t = randi([1, 100]);
x_0 = @(t) 2*t + 2;
y_0 = @(t) t + 3;
k_ad = 10^(0.3591 * log10(v_origin) + 0.0952);
k_dt = 10^(0.5441 * log10(v_opposite) - 0.0795);
r_fort = (1 + 1.34 * sqrt(k_ad^2 + k_dt^2)) * L;
r_aft = (1 + 0.67 * sqrt(k_ad^2 + k_dt^2)) * L;
r_starb = (0.2 + k_ad) * L;
r_port = (0.2 + 0.75 * k_ad) * L;
syms x y;
segmentedFunction = @(m) (m < 0) * -1 + (m >= 0) * 1;
parameter = @(x,y) power(2*x / ((1 + segmentedFunction(x)) * r_starb - (1 - segmentedFunction(x)) * r_port), 2) + ...
power(2*y / ((1 + segmentedFunction(y)) * r_fort - (1 - segmentedFunction(y)) * r_aft), 2);
fcontour(parameter, [-400,400,-400, 400], 'LineWidth', 2);
Warning: Function behaves unexpectedly on array inputs. To improve performance, properly vectorize your function to return an output with the same size and shape as the input arguments.

Accepted Answer

John D'Errico
John D'Errico on 6 Dec 2023
Edited: John D'Errico on 6 Dec 2023
As an aside, since @Cris LaPierre has already suggested the use of LevelList to do the job. And that arguably is the correct answer, at least A correct answer. But there are always multiple ways to solve a problem. In that vein, we could have used fimplicit instead.
% no need to use clc and clear in an answers script
L = 100;
v_origin = 10;
v_opposite = 20;
% There is NO need to predefine t as a random number,
% since t is just the argument to a function handle. As importantly, this
% problem does not even use x_0 and y_0. So I can just comment them out and
% still make this work. I assume you use those functions later on though.
% x_0 = @(t) 2*t + 2;
% y_0 = @(t) t + 3;
k_ad = 10^(0.3591 * log10(v_origin) + 0.0952);
k_dt = 10^(0.5441 * log10(v_opposite) - 0.0795);
r_fort = (1 + 1.34 * sqrt(k_ad^2 + k_dt^2)) * L;
r_aft = (1 + 0.67 * sqrt(k_ad^2 + k_dt^2)) * L;
r_starb = (0.2 + k_ad) * L;
r_port = (0.2 + 0.75 * k_ad) * L;
syms x y;
segmentedFunction = @(m) (m < 0) * -1 + (m >= 0) * 1;
parameter = @(x,y) power(2*x ./ ((1 + segmentedFunction(x)) .* r_starb - (1 - segmentedFunction(x)) .* r_port), 2) + ...
power(2*y ./ ((1 + segmentedFunction(y)) .* r_fort - (1 - segmentedFunction(y)) .* r_aft), 2);
So parameter is an implicit function of x and y. A contour plot fixs that function at some value or set of values, and then plot the set of points in the (x,y) plane where the function takes on exactly those given values. That is what LevelList does for us. But again, parameter defines an implicit function of x and y. And that means we can use fimplicit to do the work, to identify exactly that desired contour.
fimplicit(@(x,y) parameter(x,y) - 0.5,[-400,400,-400,400],'linewidth',2)
So we can achieve exactly the same thing as countour does, merely by recognizing what a contour is, and how it relates to an implicit function.
Contour is still a good tool in general to use though, since by showing a full set of contours, we get a much better feeling for the general shape of the surface.
  2 Comments
Cris LaPierre
Cris LaPierre on 6 Dec 2023
Very true. Another way...
L = 100;
v_origin = 10;
v_opposite = 20;
t = randi([1, 100]);
x_0 = @(t) 2*t + 2;
y_0 = @(t) t + 3;
k_ad = 10^(0.3591 * log10(v_origin) + 0.0952);
k_dt = 10^(0.5441 * log10(v_opposite) - 0.0795);
r_fort = (1 + 1.34 * sqrt(k_ad^2 + k_dt^2)) * L;
r_aft = (1 + 0.67 * sqrt(k_ad^2 + k_dt^2)) * L;
r_starb = (0.2 + k_ad) * L;
r_port = (0.2 + 0.75 * k_ad) * L;
syms x y;
segmentedFunction = @(m) (m < 0) * -1 + (m >= 0) * 1;
parameter = @(x,y) power(2*x ./ ((1 + segmentedFunction(x)) .* r_starb - (1 - segmentedFunction(x)) .* r_port), 2) + ...
power(2*y ./ ((1 + segmentedFunction(y)) .* r_fort - (1 - segmentedFunction(y)) .* r_aft), 2);
[X,Y] = meshgrid(-400:10:400);
C = parameter(X,Y);
contour(C, [0.5 0.5], 'LineWidth', 2)
子青
子青 on 8 Dec 2023
thank you sir, really helpful

Sign in to comment.

More Answers (1)

Cris LaPierre
Cris LaPierre on 6 Dec 2023
Use the LevelList Name-Value pair to specify what level(s) you want contours drawn at.
L = 100;
v_origin = 10;
v_opposite = 20;
t = randi([1, 100]);
x_0 = @(t) 2*t + 2;
y_0 = @(t) t + 3;
k_ad = 10^(0.3591 * log10(v_origin) + 0.0952);
k_dt = 10^(0.5441 * log10(v_opposite) - 0.0795);
r_fort = (1 + 1.34 * sqrt(k_ad^2 + k_dt^2)) * L;
r_aft = (1 + 0.67 * sqrt(k_ad^2 + k_dt^2)) * L;
r_starb = (0.2 + k_ad) * L;
r_port = (0.2 + 0.75 * k_ad) * L;
syms x y;
segmentedFunction = @(m) (m < 0) * -1 + (m >= 0) * 1;
parameter = @(x,y) power(2*x ./ ((1 + segmentedFunction(x)) .* r_starb - (1 - segmentedFunction(x)) .* r_port), 2) + ...
power(2*y ./ ((1 + segmentedFunction(y)) .* r_fort - (1 - segmentedFunction(y)) .* r_aft), 2);
fcontour(parameter, [-400,400,-400, 400], 'LineWidth', 2,'LevelList',0.5);
  3 Comments
Cris LaPierre
Cris LaPierre on 6 Dec 2023
You'd have to know something about your data if you are only wanting to plot a single level. You can gain that insight either visually by plotting, or artificially through the use of min and max functions. Once you have that insight, it's just a matter of picking a point to highlight.
You could create a contour plot first and capture the output, then replot using the output you captured to set the level.
L = 100;
v_origin = 10;
v_opposite = 20;
t = randi([1, 100]);
x_0 = @(t) 2*t + 2;
y_0 = @(t) t + 3;
k_ad = 10^(0.3591 * log10(v_origin) + 0.0952);
k_dt = 10^(0.5441 * log10(v_opposite) - 0.0795);
r_fort = (1 + 1.34 * sqrt(k_ad^2 + k_dt^2)) * L;
r_aft = (1 + 0.67 * sqrt(k_ad^2 + k_dt^2)) * L;
r_starb = (0.2 + k_ad) * L;
r_port = (0.2 + 0.75 * k_ad) * L;
syms x y;
segmentedFunction = @(m) (m < 0) * -1 + (m >= 0) * 1;
parameter = @(x,y) power(2*x ./ ((1 + segmentedFunction(x)) .* r_starb - (1 - segmentedFunction(x)) .* r_port), 2) + ...
power(2*y ./ ((1 + segmentedFunction(y)) .* r_fort - (1 - segmentedFunction(y)) .* r_aft), 2);
fc = fcontour(parameter, [-400,400,-400, 400], 'LineWidth', 2);
LL = fc.LevelList
LL = 1×7
0.5000 1.0000 1.5000 2.0000 2.5000 3.0000 3.5000
figure
fc = fcontour(parameter, [-400,400,-400, 400], 'LineWidth', 2,'LevelList',LL(1));
Dyuman Joshi
Dyuman Joshi on 6 Dec 2023
Edited: Dyuman Joshi on 8 Dec 2023
Thank you for your response @Cris LaPierre, though I feel like that is a vague description.
Could you specify which/what method you followed to get the value 0.5?

Sign in to comment.

Categories

Find more on Animation in Help Center and File Exchange

Community Treasure Hunt

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

Start Hunting!