why do I get an error "First input argument must be a function handle"

Hello,

I want to calculate the integral of the function F_z, like this:

nUn = 2000;

nU2n = 1500;

nU3n = 200;

nU4n = 20;

nVn = nUn;

nV2n = nU2n;

nV3n = nU3n;

nV4n = nU4n;

dt = 0.001 ; % Timestep

t = 0:dt:1; % Time interval of the vertical excitation

Om = 405 * 2 * pi / 60; % Frequency of exitation [rpm]

phi = Om * t;

n = 3; % Number of blades

F_z = -(nUn*cos(n*phi) + nU2n*cos(2*n*phi) + nU3n*cos(3*n*phi) + ...

+ nU4n*cos(4*n*phi) + nVn*sin(n*phi) + nV2n*sin(2*n*phi) + ...

+ nV3n*sin(3*n*phi) + nV4n*sin(4*n*phi));

% Plot of the evolution fo the vertical excitation

plot(t, F_z)

% Computing and plot of PSD

for f = 0:0.1:100

F_z_sqMean = f * integral((F_z).^2, 0, 1/f);

end

Answer 1

[Star Strider]

I have no idea what the code does.

See if this paproximates the desired result —

nUn = 2000; 
nU2n = 1500; 
nU3n = 200; 
nU4n = 20; 
nVn = nUn; 
nV2n = nU2n; 
nV3n = nU3n; 
nV4n = nU4n; 
dt = 0.001 ; % Timestep
t = 0:dt:1;  % Time interval of the vertical excitation
Om = 405 * 2 * pi / 60; % Frequency of exitation [rpm]
phi = @(t) Om * t;                                                      % <— Is This Part Of The Desired Outcome?
n = 3; % Number of blades 
F_z = @(t) -(nUn*cos(n*phi(t)) + nU2n*cos(2*n*phi(t)) + nU3n*cos(3*n*phi(t)) + ...
    + nU4n*cos(4*n*phi(t)) + nVn*sin(n*phi(t)) + nV2n*sin(2*n*phi(t)) + ...
    + nV3n*sin(3*n*phi(t)) + nV4n*sin(4*n*phi(t))); 

% Plot of the evolution fo the vertical excitation 
figure
plot(t, F_z(t)) 
grid
xlabel('t')
ylabel('F_z(t)')

% Computing and plot of PSD 
f = 0:0.1:100;
for k = 1:numel(f)
    F_z_sqMean(k,:) = f * integral(@(t)F_z(t).^2, 0, 1/f(k));   
end

Warning: Reached the limit on the maximum number of intervals in use. Approximate bound on error is 1.1e+33. The integral may not exist, or it may be difficult to approximate numerically to the requested accuracy.

figure
surfc(F_z_sqMean, 'EdgeColor','none')
grid on
set(gca, 'XScale','log', 'YScale','log', 'ZScale','log')
view(115,30)
xlabel('X')
ylabel('Y')

Experiment!

.