added mass and damping in VIV post processing

codes for added mass and damping in linear analysis

- the definition of added mass and damping, refer to Vikestad, 2000, jfs

addDamping.m

%% Description: for added mass and daming calculation based on liner analysis, math derivation, see k. vikestad, 2000, jfs
%%% input data: data2.txt; time|displacement|force|velocity
%%% output data: 
%%% author: Guofei, HIT
%% Some Parameters for Constant Calculations
clear;clc;
Dia=0.1;  % diameter
Len=6*Dia;  % length
mass =0.0138;  
rho=1.225; 
deltaT = 0.005; % time step
fs = 200;   
IBEGIN=4000; 
IEND = 12623; 
ksi = 0.0054;   % damping ratio
fn=1.4607;   % natural frequency
%% Data Loading 
Files = dir(fullfile('/home/chengpeng/Documents/HIT/VIV/Data/added_mass/','*.txt')); % file directory
data = load(Files.name);  %load data
time = data(:,1);  % extract first column data, i.e. time vector
dis = data(:,2);  % displacement
force = data(:,3); 
vel = data(:,4);
f = fftAnalysis(force,fs);  % load function: fftAnalysis
% use function: 'fft_filter.m'
force_filter =fft_filter(force,0,5,0.005); 
%figure;plot(time,force, 'r-');hold on; plot(time, force_filter, 'b-'); %Needed to be Corrected
dis_fluc = dis - mean(dis); 
vel_fluc = vel - mean(vel); 
force_fluc = force_filter - mean(force_filter);
f = fftAnalysis(force_fluc,fs); % call 'fftanalysis' function
acc_fluc = diff(vel_fluc)/deltaT; acc_fluc(length(vel),1) = acc_fluc(length(vel)-1,1);
acc_filter =fft_filter(acc_fluc,0,5,0.005); acc_filter =acc_filter -mean(acc_filter);
%figure;plot(time,acc_fluc, 'r-');hold on; plot(time, acc_filter, 'b-'); %Needed to be Corrected
%figure; plot(time, force_fluc, 'r-');hold on; plot(time, acc_filter, 'b-'); %Needed to be Corrected
%figure; subplot(1,2,1);plot(time, force_fluc, 'r-');subplot(1,2,2); plot(time, acc_filter, 'b-'); %Needed to be Corrected
FA = -force_fluc.*acc_filter;
FV = -force_fluc.*vel_fluc; FV = -force_fluc.*vel_fluc; FV = fft_filter(FV,0,5,0.005); 
%% Abstract Data to calculate the add-mass force and add-damping force
FAC = FA(IBEGIN:IEND);FVC = FV(IBEGIN:IEND); TIMEC = time(IBEGIN:IEND); DISC = dis_fluc(IBEGIN:IEND);
peaksFAC=findpeaks(FAC);peaksFVC=findpeaks(FVC);
IndMin=find(diff(sign(diff(FAC)))>0)+1;  IndMinFVC = find(diff(sign(diff(FVC)))>0)+1;   
IndMax=find(diff(sign(diff(FAC)))<0)+1;  IndMaxFVC = find(diff(sign(diff(FVC)))<0)+1; 
%figure;plot(FAC(1:end)); hold on; plot(IndMin,FAC(IndMin),'r^'); hold on; plot(IndMax,FAC(IndMax),'k*'); 
%figure;plot(FVC(1:end)); hold on; plot(IndMinFVC,FVC(IndMinFVC),'r^'); hold on; plot(IndMaxFVC,FVC(IndMaxFVC),'k*'); 

for i =1:length(IndMax)-1;
    aveFAC(i,1) = mean(FAC(IndMax(i):IndMax(i+1)));
    TIMEAC(i,1) = mean(TIMEC(IndMax(i):IndMax(i+1)));
    Amp = max(abs(DISC(IndMax(i):IndMax(i+1))));
    MA(i,1) = 2.*aveFAC(i,1)/(2.*pi*f(1,1))^4/Amp/Amp;
end
for i =1:length(IndMaxFVC)-1;
    aveFVC(i,1) = mean(FVC(IndMaxFVC(i):IndMaxFVC(i+1)));
    TIMEVC(i,1) = mean(TIMEC(IndMaxFVC(i):IndMaxFVC(i+1))); 
    Amp = max(abs(DISC(IndMaxFVC(i):IndMaxFVC(i+1)))); 
    CA(i,1) = 2.*aveFVC(i,1)/(2.*pi*f(1,1))^2/Amp/Amp;
end
MA_star = MA./(rho*Dia^2/4*Len);
CA_star = CA./(4.*pi*fn*mass*ksi);
figure;subplot(1,2,1);plot(TIMEAC,MA_star,'ro-');subplot(1,2,2);plot(TIMEVC,CA_star,'bs-');
Matrix_data1=zeros(length(TIMEAC),2);
Matrix_data1(:,1)=TIMEAC;
Matrix_data1(:,2)=MA_star;
Matrix_data2=zeros(length(TIMEVC),2);
Matrix_data2(:,1)=TIMEVC;
Matrix_data2(:,2)=CA_star;
save('dat1', 'Matrix_data1');
save('dat2', 'Matrix_data2');
mean(MA_star)
mean(CA_star)

fftAnalyisis.m

% FFT: Single-Sided Amplitude Spectrum of X(t)
function fd = fftAnalysis(u,fs)
L=length(u);  
Y = fft(u);
P2 = abs(Y/L); 
P1 = P2(1:L/2+1); 
P1(2:end-1) = 2*P1(2:end-1); 
f = fs*(0:(L/2))/L;
fp =[f', P1];
% fd = sortrows(fp,-2);
fd = sortrows(fp,2,'descend');
% plot(f,P1) 
% title('Single-Sided Amplitude Spectrum of X(t)')
% xlabel('f (Hz)')
% ylabel('|P1(f)|')

 fft_filter.m

function filter_data=fft_filter(input_data,low_f,high_f,delt)
%FFT_FILTER  Discrete Fourier Transform Filter
%    FILTER_DATA = FFT_FILTER (INPUT_DATA,LOW_F,HIGH_F,DELT)
%    The the vector INPUT_DATA (length N) is filtered base on the discrete
%    Fourier transform.
%
%    INPUT_DATA is the time series of length N.
%    LOW_F is the lowest frequence which can pass the filtering,
%          if  LOW_F = 0, it is the low-pass filtering.
%    HIGH_F is the highest frequence which can pass the filtering,
%           if  HIGH_F = 999, it is the high-pass filtering.
%    DELT is amount of time between each INPUT_DATA value, 
%         i.e. the sampling time.
%    The data output is also the vector of length N.
%
%----------------------------------------------------------------------------
%   Copyright (C) 2006_2009, Song Dehai
%   Ocean University of China, Program in Physical Oceanography.
%
%--------------------------------------------------------------------------
%--
if nargin~=4
     error('input error')
     return
end
%
[row_data,column_data]=size(input_data);
N=max(row_data,column_data);
if low_f==0
    fj=0;
else
    fj=low_f;
end
if high_f==999
    fk=1/delt;
else
    fk=high_f;
end
jj=max(fix(fj*N*delt),1);
kk=min(fix(fk*N*delt),N);

%Discrete Fourier transform
bf=fft(input_data-mean(input_data),N);

wm=zeros(row_data,column_data);
wm(jj:kk)=2*bf(jj:kk);

% wm([1:jj kk:fix(400*delt*N)])=2*bf([1:jj kk:fix(400*delt*N)]);
%Inverse discrete Fourier transform.
fwm=ifft(wm,N);

filter_data=mean(input_data)+real(fwm);

end

test input data

time, displacement, force, velocity

8.46	1.85E-4	1.2E-4	-0.00158
8.465	1.77E-4	1.12E-4	-0.00161
8.47	1.69E-4	1.04E-4	-0.00164
8.475	1.61E-4	9.6E-5	-0.00167
8.48	1.53E-4	8.8E-5	-0.0017
8.485	1.44E-4	7.9E-5	-0.00173
8.49	1.36E-4	7.1E-5	-0.00176
8.495	1.27E-4	6.2E-5	-0.00179
8.5	1.18E-4	5.4E-5	-0.00182
8.505	1.09E-4	4.5E-5	-0.00185
8.51	1E-4	3.6E-5	-0.00187
8.515	9.1E-5	2.7E-5	-0.0019
8.52	8.1E-5	1.8E-5	-0.00192
8.525	7.2E-5	1E-5	-0.00194
8.53	6.2E-5	1E-6	-0.00197
8.535	5.2E-5	-8E-6	-0.00199
8.54	4.2E-5	-1.7E-5	-0.00201
8.545	3.2E-5	-2.6E-5	-0.00203
8.55	2.2E-5	-3.5E-5	-0.00205
8.555	1.2E-5	-4.4E-5	-0.00207
8.56	2E-6	-5.4E-5	-0.00209
8.565	-9E-6	-6.3E-5	-0.0021
8.57	-1.9E-5	-7.2E-5	-0.00212
8.575	-3E-5	-8.1E-5	-0.00213
8.58	-4E-5	-9E-5	-0.00214
8.585	-5.1E-5	-9.9E-5	-0.00215
8.59	-6.2E-5	-1.08E-4	-0.00216
8.595	-7.2E-5	-1.17E-4	-0.00217
8.6	-8.3E-5	-1.26E-4	-0.00218
8.605	-9.4E-5	-1.35E-4	-0.00218
8.61	-1.05E-4	-1.44E-4	-0.00219
8.615	-1.16E-4	-1.53E-4	-0.00219
8.62	-1.27E-4	-1.62E-4	-0.00219
8.625	-1.38E-4	-1.71E-4	-0.0022
8.63	-1.49E-4	-1.8E-4	-0.00219
8.635	-1.6E-4	-1.89E-4	-0.00219
8.64	-1.71E-4	-1.98E-4	-0.00219
8.645	-1.82E-4	-2.07E-4	-0.00218
8.65	-1.93E-4	-2.15E-4	-0.00218
8.655	-2.04E-4	-2.24E-4	-0.00217
8.66	-2.14E-4	-2.33E-4	-0.00216
8.665	-2.25E-4	-2.41E-4	-0.00215
8.67	-2.36E-4	-2.5E-4	-0.00214
8.675	-2.47E-4	-2.58E-4	-0.00212
8.68	-2.57E-4	-2.66E-4	-0.00211
8.685	-2.68E-4	-2.75E-4	-0.00209
8.69	-2.79E-4	-2.83E-4	-0.00207
8.695	-2.89E-4	-2.91E-4	-0.00205
8.7	-2.99E-4	-2.99E-4	-0.00203
8.705	-3.09E-4	-3.07E-4	-0.00201
8.71	-3.2E-4	-3.15E-4	-0.00199
8.715	-3.3E-4	-3.23E-4	-0.00196
8.72	-3.39E-4	-3.31E-4	-0.00194
8.725	-3.49E-4	-3.38E-4	-0.00191
8.73	-3.59E-4	-3.46E-4	-0.00188
8.735	-3.68E-4	-3.54E-4	-0.00185
8.74	-3.78E-4	-3.61E-4	-0.00182
8.745	-3.87E-4	-3.68E-4	-0.00179
8.75	-3.96E-4	-3.76E-4	-0.00175
8.755	-4.05E-4	-3.83E-4	-0.00172
8.76	-4.13E-4	-3.9E-4	-0.00168
8.765	-4.22E-4	-3.97E-4	-0.00165
8.77	-4.3E-4	-4.04E-4	-0.00161
8.775	-4.38E-4	-4.11E-4	-0.00157
8.78	-4.46E-4	-4.18E-4	-0.00153
8.785	-4.54E-4	-4.24E-4	-0.00149
8.79	-4.62E-4	-4.31E-4	-0.00145
8.795	-4.69E-4	-4.38E-4	-0.00141
8.8	-4.76E-4	-4.44E-4	-0.00137
8.805	-4.83E-4	-4.5E-4	-0.00133
8.81	-4.9E-4	-4.57E-4	-0.00128
8.815	-4.96E-4	-4.63E-4	-0.00124
8.82	-5.03E-4	-4.69E-4	-0.0012
8.825	-5.09E-4	-4.75E-4	-0.00115
8.83	-5.15E-4	-4.81E-4	-0.00111
8.835	-5.2E-4	-4.86E-4	-0.00106
8.84	-5.26E-4	-4.92E-4	-0.00102
8.845	-5.31E-4	-4.98E-4	-9.74E-4
8.85	-5.36E-4	-5.03E-4	-9.29E-4
8.855	-5.41E-4	-5.08E-4	-8.84E-4
8.86	-5.45E-4	-5.13E-4	-8.39E-4
8.865	-5.5E-4	-5.18E-4	-7.94E-4
8.87	-5.54E-4	-5.23E-4	-7.5E-4
8.875	-5.57E-4	-5.28E-4	-7.05E-4
8.88	-5.61E-4	-5.33E-4	-6.61E-4
8.885	-5.65E-4	-5.38E-4	-6.17E-4
8.89	-5.68E-4	-5.42E-4	-5.73E-4
8.895	-5.71E-4	-5.46E-4	-5.3E-4
8.9	-5.73E-4	-5.51E-4	-4.87E-4
8.905	-5.76E-4	-5.55E-4	-4.45E-4
8.91	-5.78E-4	-5.59E-4	-4.03E-4
8.915	-5.8E-4	-5.62E-4	-3.61E-4
8.92	-5.82E-4	-5.66E-4	-3.21E-4
8.925	-5.84E-4	-5.7E-4	-2.8E-4
8.93	-5.86E-4	-5.73E-4	-2.41E-4
8.935	-5.87E-4	-5.76E-4	-2.02E-4
8.94	-5.88E-4	-5.79E-4	-1.64E-4
8.945	-5.89E-4	-5.82E-4	-1.27E-4
8.95	-5.9E-4	-5.85E-4	-9.1E-5
8.955	-5.9E-4	-5.88E-4	-5.5E-5
8.96	-5.91E-4	-5.91E-4	-2E-5
8.965	-5.91E-4	-5.93E-4	1.4E-5
8.97	-5.91E-4	-5.95E-4	4.7E-5
8.975	-5.91E-4	-5.98E-4	7.9E-5
8.98	-5.9E-4	-6E-4	1.1E-4
8.985	-5.9E-4	-6.02E-4	1.4E-4
8.99	-5.89E-4	-6.03E-4	1.7E-4
8.995	-5.88E-4	-6.05E-4	1.98E-4
9	-5.87E-4	-6.07E-4	2.25E-4
9.005	-5.86E-4	-6.08E-4	2.52E-4
9.01	-5.85E-4	-6.09E-4	2.77E-4
9.015	-5.84E-4	-6.1E-4	3.01E-4
9.02	-5.82E-4	-6.11E-4	3.25E-4
9.025	-5.81E-4	-6.12E-4	3.47E-4
9.03	-5.79E-4	-6.13E-4	3.68E-4
9.035	-5.77E-4	-6.14E-4	3.88E-4
9.04	-5.76E-4	-6.14E-4	4.08E-4
9.045	-5.74E-4	-6.15E-4	4.26E-4
9.05	-5.71E-4	-6.15E-4	4.43E-4
9.055	-5.69E-4	-6.15E-4	4.59E-4
9.06	-5.67E-4	-6.15E-4	4.75E-4
9.065	-5.65E-4	-6.15E-4	4.89E-4
9.07	-5.62E-4	-6.15E-4	5.03E-4
9.075	-5.6E-4	-6.14E-4	5.15E-4
9.08	-5.57E-4	-6.14E-4	5.27E-4
9.085	-5.55E-4	-6.13E-4	5.37E-4
9.09	-5.52E-4	-6.12E-4	5.47E-4
9.095	-5.49E-4	-6.12E-4	5.56E-4
9.1	-5.47E-4	-6.11E-4	5.64E-4
9.105	-5.44E-4	-6.1E-4	5.71E-4
9.11	-5.41E-4	-6.08E-4	5.78E-4
9.115	-5.38E-4	-6.07E-4	5.83E-4
9.12	-5.35E-4	-6.06E-4	5.88E-4
9.125	-5.32E-4	-6.04E-4	5.93E-4
9.13	-5.29E-4	-6.03E-4	5.96E-4
9.135	-5.26E-4	-6.01E-4	5.99E-4
9.14	-5.23E-4	-5.99E-4	6.01E-4
9.145	-5.2E-4	-5.98E-4	6.03E-4
9.15	-5.17E-4	-5.96E-4	6.04E-4
9.155	-5.14E-4	-5.94E-4	6.05E-4
9.16	-5.11E-4	-5.92E-4	6.05E-4
9.165	-5.08E-4	-5.89E-4	6.04E-4
9.17	-5.05E-4	-5.87E-4	6.03E-4
9.175	-5.02E-4	-5.85E-4	6.02E-4
9.18	-4.99E-4	-5.82E-4	6E-4
9.185	-4.96E-4	-5.8E-4	5.98E-4
9.19	-4.93E-4	-5.77E-4	5.96E-4
9.195	-4.9E-4	-5.74E-4	5.93E-4
9.2	-4.87E-4	-5.72E-4	5.91E-4
9.205	-4.84E-4	-5.69E-4	5.87E-4
9.21	-4.81E-4	-5.66E-4	5.84E-4
9.215	-4.78E-4	-5.63E-4	5.81E-4
9.22	-4.75E-4	-5.6E-4	5.77E-4
9.225	-4.73E-4	-5.57E-4	5.73E-4
9.23	-4.7E-4	-5.54E-4	5.7E-4
9.235	-4.67E-4	-5.51E-4	5.66E-4
9.24	-4.64E-4	-5.48E-4	5.62E-4
9.245	-4.61E-4	-5.44E-4	5.58E-4
9.25	-4.58E-4	-5.41E-4	5.54E-4
9.255	-4.56E-4	-5.37E-4	5.5E-4
9.26	-4.53E-4	-5.34E-4	5.47E-4
9.265	-4.5E-4	-5.31E-4	5.43E-4
9.27	-4.47E-4	-5.27E-4	5.4E-4
9.275	-4.45E-4	-5.23E-4	5.37E-4
9.28	-4.42E-4	-5.2E-4	5.34E-4
9.285	-4.39E-4	-5.16E-4	5.31E-4
9.29	-4.37E-4	-5.12E-4	5.28E-4
9.295	-4.34E-4	-5.08E-4	5.26E-4
9.3	-4.31E-4	-5.05E-4	5.24E-4