How to change non scalar variable to scalar variable?

Good day, im tryin to solve a problem in matlab following steps which were given to us. I am to find the integral of an ode with the lower limit being 0 and the upper limit being x, with x having multiple values. The solution being s or (distance) which will also have multiple values. Im trying to correct the script but I am getting an error because the x limit is not scalar. This is step 3 of 8 and im ultimately trying to complete up to step 8.
This is what I have so far thanks to the help I received:
% k = air resistance coefficient
% 𝜃 (theta) = initial angle of elevation of the ball when it was kicked
% x and y are respectively the horizontal and vertical spatial coordinates
% t is the time measured from the moment the ball is kicked
% |v| = mag_v = is the magnitude of the velocity (speed, changes with time)
global g
m = 0.45 % in kg
g = 9.81 % gravitational acceleration in m/s^2
v0 = 108 % in km/h
normv = sqrt((v0*cos(theta).^2) + (v0*sin(theta).^2))
t = 0: 0.1: 2;
theta = 30*pi/180; %initial guess for theta
k = 3 %initial guess for air resistance in N
ICs = [0 0 v0*cos(theta) v0*sin(theta)]
% v = 10
% dx/dt = x1
% dx1 / dt = x2
% mx2 = (-k)|v|(x1)
% dy/dt = y1
% dy1 / dt = y2
% my2 = (-mg)(-k)|v|(y1)
% mx2 = @(t,x)[x(1); x(2); -k*mag_v*x(1)]
% my2 = @(t,y)[y(1); y(2); -m*g -k*mag_v*y(2)]
% y(3) = dx/ dt, y(4) = dy/ dt
[T,Y] = ode45(@(t,y)fun(t,y,g,k,m),t,ICs);
%plot(T,Y(:,1:2))
Y(:,1) % X(t)
Y(:,2) % Y(t)
x = [Y(:, 1)]'
y = [Y(:, 2)]'
interp1(Y(:,1), Y(:,2), 0.7)
f = @ (x) sqrt( 1 + (y(4) ./ y(3)).^2 )
quad (f, 0, x)
%s =
function dy = fun(t,y,g,k,m)
mag_v = sqrt(y(3)^2+y(4)^2);
dy = zeros(4,1);
dy(1) = y(3);
dy(2) = y(4);
dy(3) = -k*mag_v*y(3)/m;
dy(4) = (-m*g-k*mag_v*y(4))/m;
end

12 Comments

interp1(Y(:,1), Y(:,2), 0.7)
That calculates something and displays the result, but does not store the result.
f = @ (x) sqrt( 1 + (y(4) ./ y(3)).^2 )
That anonymous function ignores its inputs and returns a constant scalar -- not something the same size as x
I used cumtrapz but im getting an error:
x = [Y(:, 1)]'
y = [Y(:, 2)]'
x1 = [0:0.05:0.7528]
f= sqrt(1 + (y(4) ./ y(3)).^2)
interp1(Y(:,1), Y(:,2), 0.7)
s = cumtrapz(x, f)
Any idea on how this could be resolved?
Q = cumtrapz(X,Y) integrates Y with respect to the coordinates or scalar spacing specified by X.
Your f is scalar. Your x is not. So you are calling cumtrapz(x, f) with a vector x and a scalar f. You do not have a scalar X and you do not have an X that has the same number of entries as Y (or the same number of columns in the 2D case). That is an error.
Why are you attempting to take the cumulative sum of a constant value?
Maybe your f should not be a scalar?
Okay, I tried changing f to a vector but I still get an error.
Isn't s given by
ds/dt = sqrt((dx/dt)^2 + (dy/dt)^2)
instead of
ds/dt = sqrt(1+((dy/dt)/(dx/dt))^2)
?
mag_v or |v| (magnitude of v) is given as: sqrt((v0*cos(theta).^2) + (v0*sin(theta).^2))
But the formula for the distance using x and y is: s = (integral) sqrt(1+ ((dy/dt)/(dx/dt))^2) dx
sorry for the confusion. The first is the magnitude of v while the second is the distances based on the values of x and y.
If you have
x = [Y(:, 1)]'
y = [Y(:, 2)]'
then dx is diff(x) and dy is diff(y) -- but I recommed instead using gradient
Hi I have this code and following along with the problem sheet attached in the first comment, I would like to create a while loop for the least squares method at the end. Any way this could be done using a while loop?
% compute_drag_script
D = [ 0; 0.5; 1.1; 1.7; 2.3; 2.9; 3.5; 4.0; 4.6; 5.2;...
5.7; 6.3; 6.8; 7.4; 7.9; 8.4; 9.0; 9.5; 10.0; 10.6;...
11.1; 11.6; 12.1; 12.6; 13.1; 13.6; 14.1; 14.6; 15.1; 15.6;...
16.1; 16.5; 17.0; 17.5; 18.0; 18.4; 18.9; 19.4];
V = [108; 107; 106; 105; 104; 104; 103; 102; 101; 100;...
99; 99; 98; 97; 97; 96; 95; 94; 94; 93;...
92; 92; 91; 91; 90; 89; 89; 88; 88; 87;...
87; 86; 86; 85; 85; 84; 83; 82];
interp_speeds = compute_drag(D,V,0.0095,30);
function interp_speeds = compute_drag(...
true_distances,true_speeds, k,theta_degrees)
arguments
true_distances(:,1) {mustBeNonnegative}
true_speeds(:,1) {mustBeNonnegative}
k (1,1) {mustBeNonnegative}
%drag_coeff
theta_degrees(1,1) {mustBeInRange(theta_degrees,0,180)}
end
%% Input validation
if length(true_distances) ~= length(true_speeds)
error('Distance and speed vectors must be the same length.');
end
if any(diff(true_distances) < 0)
error('Distance should never go down!');
end
%% Constants and parameters
g = 9.81; % meters per second-squared
m = 0.45; % kilograms
v0_kph = 108; % kilometers per hour, will convert to meters per second
k = 0.00000001;
ro = 1.293 % kg/(m)^3
d = 0.7/pi % m
%% Conversions
v0 = v0_kph * 1000 / 3600; % kph to meters per second
% We could also convert theta_guess to radians, but there's no reason we
% have to, since Matlab has degree variants of the trig functions.
%% Initial Conditions
x0 = 0;
y0 = 0;
dxdt0 = v0*cosd(theta_degrees);
dydt0 = v0*sind(theta_degrees);
init_vals = [ x0; y0; dxdt0; dydt0];
%% Step 2
% setup for diffeq that will be fed to ode45:
% val(1) x: x' (i.e. val(3))
% val(2) y: y' (i.e. val(4))
% val(3) x': (-k|v|/m) * x'
% val(4) y': -g - (k|v|/m) * y'
all_diffeq = @(t,VALS) [
VALS(3);
VALS(4);
-k * sqrt( VALS(3)^2 + VALS(4)^2 ) * VALS(3) / m;
-g - ( k*sqrt(VALS(3)^2+VALS(4)^2) * VALS(4) / m);
];
% To get an appropriate time range, I'm using the worst case scenario:
% the slowest speed (e.g. 82 kph) taking the ball the full distance
% (e.g. 19.4 meters).
minspeed = min(true_speeds) * 1000 / 3600;
max_necessary_time = (true_distances(end)) / minspeed;
tspan = [0, max_necessary_time];
[t,vals] = ode45(all_diffeq, tspan, init_vals);
x = vals(:,1);
y = vals(:,2);
dxdt = vals(:,3);
dydt = vals(:,4);
% s = sqrt(x.^2 + y.^2);
%% Step 3
s_segments = sqrt(1 + (dydt./dxdt).^2);
s_step = [0;diff(x)];
s = cumsum( s_segments .* s_step );
v = sqrt(dxdt.^2 + dydt.^2);
v_kph = v * 3600 / 1000;
%% Step 4
interp_speeds = interp1(s,v_kph,true_distances);
if any(isnan(interp_speeds))
warning('NaN encountered. You should increase upper bound on time.');
end
end
% Step 5
while V_interpolated >= 0
E = cumsum(V_interpolated - V_var).^2
end
Hi,
Yes, a least squares optimization can be implemented using a while loop.
In the following points, I explain the steps that can be leveraged to achieve the goal.
  • Setting while loop stopping condition: I suggest setting the while loop stop condition according to an error tolerance value instead of stopping the while loop at a zero value. The below provides a sample code for the same:
tol = 0.00001; % Define a tolerance value for the error to stop optimization
while (dE_dtheta > tol || dE_dk > tol)
% Implement the algorithm and variable correction here.
end
  • Implementing the least square optimization algorithm: The algorithm is already specified in the attached PDF document under the “LEAST SQUARES METHOD” heading. If you face specific issue in the implementation, feel free to post here.
  • Correcting the variable: The correction in initial guess variables (“k” and “theta” in this case) with each "while" loop iteration is an important step of optimization. You can use gradient descent in the following way to perform the same:
alpha = 0.001; % Defines the learning rate or step size for gradient descent.
% After calculating error values in current iteration, update the optimization variables as given below.
theta = theta - alpha*dE_dtheta;
k = k - alpha*dE_dk;
while V_interpolated >= 0
E = cumsum(V_interpolated - V_var).^2
end
Your while loop must update at least one of the variables that are being tested in the condition, or else you have an infinite loop.
Exception: it is valid to have a loop similar to
while true
do some calculation
if condition is met
break;
end
end
That is, updating the variable or value being tested in a while loop is not strictly necessary if you have a test and a break statement.

Sign in to comment.

Answers (0)

Categories

Find more on Programming in Help Center and File Exchange

Asked:

on 21 Apr 2023

Commented:

on 19 May 2023

Community Treasure Hunt

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

Start Hunting!