Dimensionless variables of the differential equation

Hello all,

I hope you are all doing well. I have established a 2DOF differetial equation. Basically, the equations are the 2 spring-mass-damper system. But I am wondering how to define the variable X to be a dimentionless variable, such as X/L. In my equations, the variable X should be dimensionless but when I exclude the force (F_o the dimensionless force), the results of the equation are not dimensionless which I have already checked. My codes are attached. I think you can easily run it.

I appreciate your support.

Best wishes,

Yu

clc;
clear all;

tspan = 0:0.25:200;

% Units:
% stiffness: N/m, mass:kg, damping(c_1):Ns/m
% displacement(X):m
% alpha (stiffness ratio):k_1/k_2
% gamma (dimension ratio):a/L (structure length)
% f_opt (frequency ratio): f_1/f_2

X0 = [0 0 0 0];
%parameters-------
mu = 0.02; % mass ratio
f_opt = 1/(1+mu); %frequency ratio
xi_2 = sqrt(3*mu/8*(1+mu));
Omega_1 = 0.188; % 0.03 Hz k_1^2/m_1^2 (stiffness^2/mass^2)
Omega_2 = Omega_1*f_opt; % frequency of the second DOF k_2^2/m_2^2 (stiffness^2/mass^2)
xi_1 = 0.01; % damping ratio of first DOF xi_1=c_1/2*m_1*Omega_1

%matrix--------
M = [1+mu mu;
     1     1];
C = [2*xi_1*Omega_1 0;
     0         2*xi_2*Omega_2];
K = [Omega_1^2 0;
     0         Omega_2^2];


% state space model-------------------

% M: mass matrix, K:stiffness matrix, C:damping matrix
O = zeros(2,2);
I = eye(2);
A = [O I; -inv(M)*K -inv(M)*C];
B = [O; inv(M)];
E = [I O];
D = zeros(2,2);

% solve the equations-------------------

options = odeset('RelTol',1e-10,'AbsTol',1e-10);

[t,X] = ode45(@(t,X) QZSdamper(t,A,B,X),tspan,X0,options);

x = X(:,1:4);

% F_o = 2*Omega_2.^2.*alpha.*(sqrt(1-gamma.^2)+X(:,2)).*(-1+sqrt(1/(X(:,2).^2+2.*sqrt(1-gamma^2).*X(:,2)+1)));
% 
% f_s = Omega_2.^2*(X(:,2)-2*alpha*(sqrt(1-gamma^2)+X(:,2))*(-1+sqrt(1/(X(:,2).^2+2*sqrt(1-gamma^2)*X(:,2)+1))));
% 
% f_1 = f_s(:,1);
% 
% figure,
% plot(t,f_1);


figure,
plot(t,x(:,1)); xlabel('Time/s'),ylabel('Dimensionless displacement of primary structure')
figure,
plot(t,(x(:,2)))

sys = ss(A,B,E,D);
figure,
bodeplot(sys(1,1)) % from input 1 to output 3
bp.FrequencyScale = "linear";

[mag,phase,wout] = bode(sys(1,1));

mag = squeeze(mag);
phase = squeeze(phase);
fout = wout/(2*pi);
BodeTable = table(fout,mag,phase);

function dXdt = QZSdamper(t,A,B,X)


mu = 0.025;
alpha = 1;
gamma = 2*alpha/(2*alpha+1);

f_opt = 1/(1+mu);
Omega_1 = 0.188; % 0.03 Hz
Omega_2 = Omega_1*f_opt;

w = 0.180;

F = 0.001*sin(w*t);

F_o = 2*Omega_2^2*alpha*(sqrt(1-gamma^2)+X(2))*(sqrt(1/(X(2)^2+2*sqrt(1-gamma^2)*X(2)+1))-1);

% F_o = Omega_2^2*(sqrt(1-gamma^2)+2*alpha*(X(2))*(1/((X(2)-sqrt(1-gamma^2))^2+2*sqrt(1-gamma^2)*(X(2)-sqrt(1-gamma^2))+1)^(1/2)-1))

F = [F;F_o];

dXdt = A*X+B*(F);

end

[John D'Errico]

I'm sorry, but MATLAB has no concept of dimensionless variables. A number is just a number. The units must be applied by the user, if there are any units. If you want to pose something in terms of meters, nanometers, inches, furlongs per fortnight, a unit-less number such as a Reynolds number, or a unit-less distance, in the form of a relative distance like X/L, MATLAB does not care. It is the user who must insure what they are doing makes mathematical sense, as well as physical sense.

If you claim the results of your computation are not dimensionless, this just means you have made a mistake in the equations you formulated.

[Yu]

@John D'Errico Hello, many thanks for your comments. I think I understand your explanation reagarding the MATLAB principle.

But I am still wondering whether MATLAB can solve the tuned mass damper FRF equations. (I am not sure if you are familiar with this, there fore, I attach the pic) In this equation, all the parameters are dimensionless.

Best wishes,

Yu

[John D'Errico]

Of course it "can" solve them. If you put in the correct system of equations, ones that have a valid, well-posed solution, then you will get a valid answer.

[Torsten]

@Yu

You can use units with the Symbolic Math Toolbox. But it seems this feature is only available in symbolic computations - thus not applicable in your code above.

https://uk.mathworks.com/help/symbolic/units-of-measurement-tutorial.html

[Yu]

@John D'Errico I am sorry whether I understand correctly. It means that MATLAB can solve the equations consisting of reasonable dimensionless parameters. BUT, it can not make the 'solution' dimensionless.

[Yu]

@Torsten @John D'Errico Hello, thank you for your suggestion. I am thinking that the 'F_o' is dimensionless after a series simplification. Therefore, I can regard the 'solution' of my equations is dimensionless. But if F_o is removed the solution is not dimensionless. Is this correct?

[Torsten]

We have no clue what equations you try to solve in your code, what your parameters mean and what units they have. The included picture doesn't help in this respect because the relationship between your code and the parameters and equations used in the article/book are not obvious.

[John D'Errico]

MATLAB cannot "make" a solution dimensionless, since again, MATLAB does not understand units on variables. (With the exception of computations done with the symbolic toolbox as has been pointed out, and since you are using ODE45, that does not apply. It is also why I did not mention that capability.) Only you can insure something is dimensionless, which requires that you understand the formulation of your equations, and those equations were properly rendered dimensionless when you created them. It seems the part that is lacking is your understanding here.
MATLAB is just a tool. It does what you tell it to do, nothing more and nothing less. And it you tell it to do something physically meaningless, it will still proceed and do as it is told. The result will of course, be meaningless. But that is not the fault of MATLAB, but the user.

[Yu]

@John D'Errico@Torsten I am really sorry that my question has not clear comments. For my codes, I have modified my question regarding the contents. Also I attached the dimensionless parameters in my equations. Could you please help me with this problem when you are available? I would appreciate it.

The reason why I attached the pic before is that the equation in the textbook is dimensionless as you would see. But as my equation has the nonlinear term. So it would be difficult to transfer my equation to the same form as the equation in the textbook. Thus, I kindly ask for some help.

[Sam Chak]

Hi @Yu

I ran the code, and the dimensionless displacement becomes unstable due to disturbance from the oscillating force (in QZSdamper). For your information, F_o​ can stabilize the system but is insufficient to handle the oscillating force.

What is the expected outcome?

tspan   = 0:0.25:200;

% Units:
% stiffness: N/m, mass:kg, damping(c_1):Ns/m
% displacement(X):m
% alpha (stiffness ratio):k_1/k_2
% gamma (dimension ratio):a/L (structure length)
% f_opt (frequency ratio): f_1/f_2


X0      = [0; 0; 0; 0];
% X0      = [1; 0; 0; 0];     % test if F_o can stabilize when sinusoidal force, F = 0.

% Parameters-------
mu      = 0.02; % mass ratio
f_opt   = 1/(1 + mu);         % frequency ratio
xi_2    = sqrt(3*mu/8*(1 + mu));
Omega_1 = 0.188;            % 0.03 Hz k_1^2/m_1^2 (stiffness^2/mass^2)
Omega_2 = Omega_1*f_opt;    % frequency of the second DOF k_2^2/m_2^2 (stiffness^2/mass^2)
xi_1    = 0.01;             % damping ratio of first DOF xi_1=c_1/2*m_1*Omega_1

% Matrices 
% M: mass matrix, 
M       = [1+mu, mu;
           1   ,  1];

% C: damping matrix
C       = [2*xi_1*Omega_1, 0;
           0             , 2*xi_2*Omega_2];

% K: stiffness matrix, 
K       = [Omega_1^2, 0;
           0        , Omega_2^2];

% state space model-------------------
O   = zeros(2,2);
I   = eye(2);
A   = [O,         I; 
      -inv(M)*K, -inv(M)*C]

A = 4×4

0 0 1.0000 0 0 0 0 1.0000 -0.0353 0.0007 -0.0038 0.0006 0.0353 -0.0347 0.0038 -0.0329

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

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

ev  = eig(A)    % eigenvalues, stable but easily becomes unstable when disturbed

ev = 
  -0.0099 + 0.1972i
  -0.0099 - 0.1972i
  -0.0085 + 0.1753i
  -0.0085 - 0.1753i

B   = [O; 
       inv(M)];
E   = [I O];
D   = zeros(2,2);

% solve the equations-------------------
options = odeset('RelTol', 1e-10, 'AbsTol', 1e-10);
[t, x]  = ode45(@(t, X) QZSdamper(t, A, B, X), tspan, X0, options);

% x = X(:,1:4);

% F_o = 2*Omega_2.^2.*alpha.*(sqrt(1-gamma.^2)+X(:,2)).*(-1+sqrt(1/(X(:,2).^2+2.*sqrt(1-gamma^2).*X(:,2)+1)));
% 
% f_s = Omega_2.^2*(X(:,2)-2*alpha*(sqrt(1-gamma^2)+X(:,2))*(-1+sqrt(1/(X(:,2).^2+2*sqrt(1-gamma^2)*X(:,2)+1))));
% 
% f_1 = f_s(:,1);
% 
% figure,
% plot(t,f_1);

figure,
plot(t, x(:,1)); grid on
xlabel('Time/s'),
ylabel('Dimensionless displacement of primary structure')

% figure,
% plot(t,(x(:,2)))

% sys = ss(A,B,E,D);
% figure,
% bodeplot(sys(1,1)) % from input 1 to output 3
% bp.FrequencyScale = "linear";
% 
% [mag,phase,wout] = bode(sys(1,1));
% 
% mag = squeeze(mag);
% phase = squeeze(phase);
% fout = wout/(2*pi);
% BodeTable = table(fout,mag,phase);

function dXdt = QZSdamper(t,A,B,X)
    mu      = 0.025;
    alpha   = 1;
    gamma   = 2*alpha/(2*alpha+1);
    f_opt   = 1/(1+mu);
    Omega_1 = 0.188; % 0.03 Hz
    Omega_2 = Omega_1*f_opt;
    w       = 0.180;
    
    % oscillating force
    F       = 0.001*sin(w*t);
    
    % Unknown force (comes from somewhere)
    F_o     = 2*Omega_2^2*alpha*(sqrt(1 - gamma^2) + X(2))*(sqrt(1/(X(2)^2 + 2*sqrt(1 - gamma^2)*X(2) + 1)) - 1);
    % F_o = Omega_2^2*(sqrt(1-gamma^2)+2*alpha*(X(2))*(1/((X(2)-sqrt(1-gamma^2))^2+2*sqrt(1-gamma^2)*(X(2)-sqrt(1-gamma^2))+1)^(1/2)-1))
    
    F       = [F;       % same F = F? bad coding practice from Engineering perspective
               F_o];    % ?
    
    % F       = [0;     % no force to test stability
    %            F_o];    % no force
    
    dXdt    = A*X + B*F;
end

[Torsten]

@Sam Chak

I think the OP's main question is whether her/his equations are per se correct because of the "inconsistency with F_o". But it's still unclear what she/he tries to explain.

[Yu]

@Sam Chak@Torsten Hi, thank you for your comments. I think I would explain my model in more details. I am not sure whether you have time to review it. But I suppose they could be helpful.

My model is the 2-DOF system, including the primary structure and the damper (seen below)

Therefore, the equations:

f_o is the nonlinear stiffness force from oblique spring the strcuture below. k_tmd i from the vertical spring. So if f_o is removed, it will become the tuned mass damper not QZS.

Fig.1

The above structure (Fig.1) can be expressed by (k_v vertical spring stiffness, k_o oblique spring)

Then, the first eqn is divided by m_1, and second is divided by m_tmd. The eqns are now:

Actually, the above equations are checked in symbolic tools. As the force generate by the QZS (Fig.1) is related to the dimensional parameters, therefore, I just want to use the ratios as I mentioned in the questions. I do not want to use the specific numbers at this moment.

Regarding the unstable problem, I removed the f_o and the tspan is 1000. The system are stable? Therefore, I guess my expression may be correct? Also if you see the bodeplot, there is a peak reduction I believe and the structural damping ratio is also low.

For employing F_o (pic below), it seems it not that stable. I guess the stiffness has a 'zero' point that there is a position the total stiffness is 0. For the equations of the F_o, I also verified. It seems it is correct.

The final question is that I am not sure the results of the displacment is dimensionless or not. Especially when removing F_o.

[John D'Errico]

I'm starting to see what you did wrong. Almost enough for my response here to qualify as an answer, but not really.
You divided the equations by length, thinking this makes the length a dimensionless variable. But force has units of mass*distance/time^2. (For example, a common unit of force might be represented as Kg*M/sec^2.)
It is not at all clear what you did to make force dimensionless. Dividing by length does not result in a dimensionless force.

[Torsten]

How are the vectors f_o and F(t) in your very first equation defined ? Is f_o = [0 ;f_non] ? Is F(t) = [a*sin(omega*t);0] ?

Is f_o~ = fbar ?

Is x_tmd = x_2 ?

Is x = [x_1;x_2] in the two grey graphics or are x_1 and x_2 also non-dimensionalized in the second grey graphics ?

There is a hopeless mess of different variable names in your description.

Don't you have to present your model together with your results to your professor in a transparent way ?

[Torsten]

@John D'Errico

I guess displacement and time must be non-dimensionalized, too. Maybe it's the case in the second gray graphics.

The transformation should somehow be

x~(t~) = 1/L * x(t)  

for some reference length L and end time T with t = t~ * T.

x~ and t~ are dimensionless here.

[John D'Errico]

@Torsten

It is not sufficient to non-dimensionalize only displacement and time, since the units on that mess are units of force.

I tried to find where that was done, but I could not. All that was said is repeatedly how distance was non-dimensionalized. And that is insufficient. Since the problem apparently lies with force, it is my conjecture this is where the problem lies. But that is a completely wild guess.

Again, force would have units of mass*distance/time^2.

[Yu]

@John D'Erric@Torsten Hi, thank you for your comments and explantions.

I would clarify further about my equations:

Let me assume that only F_o is here, the equation is

Then, (4) dividede by kv*L_o, it woul be (here is a relationship: gamma=a/L_0). k_v is just used for eliminating the k element to achieve dimensionless form. alpha = k_o/k_v

Then, if I use this approach in the second equation

The stiffness force in Eq. (6) consists of linear and nonlinear terms. Then, I divided m_tmd both sides

I think here is clear. Then to get to the Eq. (5), I divided L_0 on both sides

Now, I can regard x2/L_0 is my results. Are these of my thoughts correct?

I would appreciate your kindly help.

[Torsten]

@John D'Errico

Again, force would have units of mass*distance/time^2.

Yes. Dividing each equation by m and the constant L/T^2 (which results from the non-dimensionalization of displacement and time by x~(t~) = 1/L * x(t) ) makes force dimensionless.

[Yu]

@Torsten I am really sorry that I have presented a mess description of my questions. I made some screenshots from word.doc last night when I got the comments from you.

I am aware that it may waste your time to undestand the mess variables and apologize for this.

I have editted my previous reply. I also clarify here

1. F(t) is the harmonic excitation on the m_1. It should be [f(t) 0].
2. Yes, F(t) is asin(omega*t) --- the harmonic excitation
3. f_o = . This also is modified in my previous reply.
4. Yes, x_tmd = x_2
5. I just nondimensionalize the nonlinear stiffness force f_o in the second equation. But I have no clue about how to make the first equation nondimensionalized. Additionally, I am not sure I nondimensionalize the second equation successfully.
6. Yes I need to present the clear description of my question to the professor. It was my mistake.

Again, I am sorry for the inconvenience and thank you for your patient help.

Best wishes,

Yu

Answer 1

[Sam Chak]

Hi @Yu

The aim is to render the entire dynamical system dimensionless (by finding the dimensionless time derivative), rather than focusing solely on a single state. I didn't study your tuned mass damper with with quasi-zero-stiffness (QZS), so I provide a relatively simple example below:

Example:

Consider the following example, which presents a simplified set of coupled dynamic equations for a tethered space system (TSS):

where

• l is the tether length (unit: m),
• θ is the libration angle (unit: radians, but considered dimensionless),
• T is the tether tension (unit: kg·m/s²)
• ω is the orbital angular velocity (a constant) of TSS (unit: rad/s²), and
• is the mass equivalent (a constant), analogous to adding two "parallel resistors"; the mothership and the tethered probe (unit: kg).

To non-dimensionalize the dynamics of tether length, consider making the tether tension dimensionless by dividing the entire equation by a constant, . Note that is a nominal tether length chosen to satisfy a physical constraint.

If the following dimensionless quantities are defined:

• dimensionless tension,
• dimensionless length, , such that
• dimensionless time,
• dimensionless time derivative,

then the dimensionless form of the TSS dynamic equations can written as:

[Yu]

Hi Sam,

Many thanks for your answer. It is really clear about the nondinmensionlize the equations.

May I ask a question?

in your statement is the nominal tether. But may I know how to determine it?

For instance, in my equations, , (α is determined here) and . It means I need to determine the real values of these three terms that can statisfy the relationship above:

e.g.

m, m, and m. By dividing , I can also reach out the dimensionless form.

Best wishes,

Yu

[Sam Chak]

Hi @Yu

Please do not confuse the TSS equations with your TMD equations; they are two distinct systems, though both share Newton's . I reviewed the "scattered equations" (which was very time-consuming), but I did not find the nominal tether length or the orbital angular velocity in your work. You should determine your own nominal displacement and the angular velocity of the TMD system.

I suggest consolidating all equations into a single post to enhance readability. This approach will also facilitate your explanation of the modifications made to the equations and make it easier for others to review them.

[Yu]

@Sam Chak Hi Sam,

Many thanks for your suggestions. And I am sorry that the unclear desciption have taken your much time. I have posted another questions with the description of my system. I guess it is clearer.