Time vs displacement plot of a Transfer function

I am trying to plot time vs displacement plot of a transfer function. It is a simple case of spherical body in a liquid hence excitation of a spherical body is being resisted by the (elastic and viscous properties of the liquid) Hence F excitation= F inertia + F (elastic) + F viscous.

The final transfer function is of the form below where Delta Q is the relative displacement between the body and the liquid where as Vm is the excitation velocity.

I have used this transfer function to find eigen frequency for this system. Now i wants to check if the system undergoes resonance at that eigen frequency. How to do it using matlab

[Sam Chak]

@Hassan

Please type out the MATLAB code and provide the parameters.

%% parameters
R   = ...;
zeta= ...;
iota= ...;
rhom= ...;
BoK = ...;
B1u = ...;
B1K = ...;

%% coefficients in numerator
n0  = 0;
n1  = 4/3*pi*(R^3)*zeta*iota*rhom;

%% coefficients in denominator
d0  = BoK;
d1  = B1u + B1K;
d2  = ...;

%% transfer function
num = [n1 n0];
den = [d2 d1 d0];
G   = tf(num, den)

[Hassan]

%% parameters
R   = 5e-6;
zeta= 0.2308;
iota= 1.3;
rhom= 997;
rhonu = 1.3*rhom;
BoK = 8.54E-01 ;
B1u = 6.03E-07 ;
B1K = 1.42E-06 ;
B2u = 6.65E-13; 
B2k = 1.31E-13; 

%% coefficients in numerator
n0  = 0;
n1  = 4/3*pi*(R^3)*zeta*iota*rhom;

%% coefficients in denominator
d0  = BoK;
d1  = B1u + B1K;
d2  = B2u+B2k+((4/3)*3.1415*(R^3)*rhom);
    
%% transfer function
num = [n1 n0];
den = [d2 d1 d0];
G   = tf(num, den)

G =
 
              1.566e-13 s
  -----------------------------------
  1.318e-12 s^2 + 2.023e-06 s + 0.854
 
Continuous-time transfer function.

[Hassan]

I have used this, and got the transfer function

4.949e-15 s

----------------------------------

1.354e-12 s^2 + 1.85e-06 s + 0.854

[Hassan]

How to get a time vs amplitude plot if;

Delta Q is output and Vm is input, input is sinusoidal at the eigen frequency of the system

[Sam Chak]

@Hassan, click the indentation icon to toggle the MATLAB coding mode. Then click the play icon to run the code. I edited the coding format in your post and ran it, but it threw an error. There are missing info.

[Hassan]

Hi Sam,

Please check, its working now.

[Sam Chak]

Hi @Hassan

Your computed transfer function is different for the one displayed in your comment. Would you check the details?

Also, please specify the simulation time window and provide the sinusoidal input signal Vm.

Answer 1

[sai charan sampara]

Hello Hassan,

The following code might help you:

%% parameters
R   = 5e-6;
zeta= 0.2308;
iota= 1.3;
rhom= 997;
rhonu = 1.3*rhom;
BoK = 8.54E-01 ;
B1u = 6.03E-07 ;
B1K = 1.42E-06 ;
B2u = 6.65E-13; 
B2k = 1.31E-13; 

%% coefficients in numerator
n0  = 0;
n1  = 4/3*pi*(R^3)*zeta*iota*rhom;

%% coefficients in denominator
d0  = BoK;
d1  = B1u + B1K;
d2  = B2u + B2k + (4/3*pi*(R^3)*rhonu);
    
%% transfer function
num = [n1 n0];
den = [d2 d1 d0];
G   = tf(num, den)

G =
 
              1.566e-13 s
  -----------------------------------
  1.475e-12 s^2 + 2.023e-06 s + 0.854
 
Continuous-time transfer function.

This gives the expression for the transfer function "G". Once the transfer function is obtained you can use the following code to obtain the eigen frequencies:

eigen_values=eig(G);
eigen_freq=abs(imag(eigen_values))

eigen_freq = 2x1

1.0e+05 * 3.2958 3.2958

<mw-icon class=""></mw-icon>

<mw-icon class=""></mw-icon>

Once the eigen frequencies are obtained you can define the input signal as a sinusoid of the eigen frequency and pass it through the system using the "lsim" function in MATLAB. This gives the output signal and the required Amplitude vs Time plot.

You can learn more about the function "lsim" from this documentation:

https://www.mathworks.com/help/control/ref/dynamicsystem.lsim.html

[Hassan]

Thankyou Sam. There is one more issue, that values like B1U and B2U are dependent on frequency itself. is there anyway i can share the details of my program in private and get your opinion. as i cannot share here, it is the part of my ongoing project.

For the values i shared earlier, I found B1U and B2U after many trials to find resonant frequency of this system. I give frequency sweep signal from 10^3 to 10^ 7, and started to narrow down until I find almost accurate frequency which gave max amplitude. Is there any other way to find it directly?

[Hassan]

Moreover, If i apply force of U=asinw(wt). where w= eigen freq., you shown me earlier, it should vibrate at increasing amplitude or growing amplitude with the time. isnt