道路动态称重的精细逐步积分法

pumac

各位好，我是新来的。

当车辆以20km/h的速度压过动态称台时，会记录下一连串的非线性的压力信号，一般包括车辆振动和车自重。如何用精细逐步积分法来分析获得的非线性压力曲线，并得到一个大致压力的结果。

能否介绍些精细逐步积分法的内容？

风花雪月

关于动态称重的理论分析可以参考
http://www.quzhi.net/chinese/pdf/WIM.pdf

风花雪月

不知道你想用精细逐步积分法，是否是求解动力学方程？

至于其原理可以参考钟万勰的论著

[ 本帖最后由 风花雪月 于 2007-7-2 10:55 编辑 ]

风花雪月

下面是一个精细积分求解结构动力学方程的程序，你可以参考一下

1. %clear;
   
2. %A=zeros(2);
   
3. %C=A;
   
4. %D=[0.5,0;0,1];
   
5. %B=[-6,2;2,-4];
   
6. %定义质量矩阵
   
7. M = [1 0 0;0 1 0;0 0 6];
   
8. %定义刚度矩阵
   
9. K = [10 1 0; 1 10 1;0 1 10];
   
10. %定义比例阻尼
    
11. alpha = 0.01;
    
12. bita = 0.02;
    
13. G = alpha*M+bita*K;
    
14. iM = inv(M);
    
15. A = -0.5*iM*G;
    
16. B = 0.25*G*iM*G-K;
    
17. C = -0.5*G*iM;
    
18. D = iM;
    
19. %
    
20. %A=zeros(3);
    
21. %C=A;
    
22. %D=[0.5 0 0;0 1 0;0 0 2];
    
23. %B=[-6 2 2;2 -4 2; 2 2 -5];
    
24. %定义初始时刻的外力向量
    
25. %f0=[0;0;0;10];
    
26. f0=[0 0 0 10 1 20]';
    
27. f1=zeros(size(f0));
    
28. %形成哈密尔顿矩阵
    
29. H=[A,D;B,C];
    
30. %计算的结束时间
    
31. tf=20;
    
32. %下面的定义是为了验证在不同步长情况下，精细积分的精度之间的差别
    
33. step=[2,0.5,0.1]; % different step size
    
34. %这一般是精细积分计算的缺省值
    
35. N=20;
    
36. %figure;
    
37. grid;
    
38. hold on;
    
39. str=['o','x','b-'];
    
40. for jj = 1:3
    
41.    [jj tf/step(jj)]  %time,vector的维数和tf/step(jj)的值相等
    
42.    [time,vector]=jxjf1(H,f0,f1,step(jj),tf,N);
    
43.    plot(time,vector(:,1),'r-');
    
44.    plot(time,vector(:,2),'g-');
    
45.    plot(time,vector(:,3),'b-');
    
46. end
    
47. %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    
48. function [time,matrix]=jxjf1(H,f0,f1,time_step,tf,N)
    
49. %
    
50. % time:计算的时间序列,一列
    
51. % matrix:计算的结果，矩阵，按照列分为两半，前一半为节点位移的响应，前一半为节点动量的响应，
    
52. %
    
53.     %定义与哈密尔顿矩阵同阶的单位矩阵
    
54.     I=eye(size(H));
    
55.    
    
56.     %需要计算增扩矩阵的逆
    
57.     iH=inv(H);
    
58.    
    
59.     %精细积分的输入为：步长、一个正整数、增扩矩阵、外力向量初值
    
60.     %PIM begin 精细积分方法开始
    
61.     dt=time_step/2^N;  %将每一个步长分为2的整幂次方份
    
62.    
    
63.     % Ta=H*dt+(H*dt)^2*(I+(H*dt)/3+(H*dt)^2/12)/2;
    
64.     temp = H*dt;
    
65.     temp1 = temp*temp;
    
66.     Ta=temp+temp1*(I+temp/3+temp1/12)/2;
    
67.    
    
68.     %对每一份计算
    
69.     for iter=1:N  
    
70.         Ta=2*Ta+Ta*Ta;  %化简成Ta=（2+Ta）*Ta有些问题 huhu
    
71.     end
    
72.     %定义一个新变量
    
73.     T=I+Ta;
    
74.    
    
75.     %vk=[0;0;0;0];  %定义一个列向量
    
76.     vk=zeros(max(size(H)),1);  %定义一个列向量 huhu
    
77.    
    
78.     number = tf/time_step;
    
79.     %number
    
80.     matrix = [];
    
81.     for iter=1:number;
    
82.       
    
83.         matrix = [matrix vk];
    
84.       
    
85.         %iter
    
86.         t(iter)=time_step*(iter-1);
    
87.         %v(iter)=vk(1);
    
88.         %vk=T*(vk+iH*(f0+iH*f1))-iH*(f0+iH*f1+f1*step(jj));
    
89.         temp = iH*(f0+iH*f1);
    
90.         vk=T*(vk+temp)-temp-iH*f1*time_step;
    
91.     end
    
92.     time = t(1:number)';
    
93.     %vector = v;
    
94.     matrix = matrix';
    
95. end

复制代码