为LAOS波形重构做了一个GUI。
波形重构GUI
除了上一篇文章提到的一些问题外，还有其他的难题。最主要的就是，当原信号的数据量特别大的时候，正弦波拟合函数和MITLaos都会出错。我设计了一个判断，适当时候把信号down sample了再做正弦波拟合和MITLaos，前者基本没问题了，但后者还是不时出问题。只要正弦波拟合没问题的话，波形重构是没问题的，没必要卡在MITLaos上。所以我打算不用MITLaos，直接在重构了的波形上求导得到G_M、G_L等参数。事实上，在程序中G_k参数就是通过直接求导得到的，结果令人满意。
对真实世界的信号进行求导，需要先大幅提高信号的信噪比，否则求导运算就是一个噪音放大器。对于LAOS实验的信号，由于假设其为一系列奇次谐波的和：
。因此，可以据此直接进行波形重构。即通过FFT得到信号的谐波振幅和相位角，按照此式重构一个平滑的波型，再拿去做微分。原理很简单，但是在MATLAB里整了半天，强度是对的，但相位角永远是错的。研究了一晚上，才发现FFT得到的相位角总是少了90°，也不知道什么原因。总之加上这90°之后就完美了。为了减少日后摸索，把代码贴在这里：
clear all;close all;clc;
%============ Generate the signal for reconstruction ============%
fs=57; % sampling frequency
fi0=10/3; % input frequency
t=0:1/fs:123456/fs; % the length of signal is 12345.
t=t';
delta=[1 3 2 1];
x1=sin(2*pi*fi0*t); % the input signal is a sinusoidal wave;
x2=12*sin(2*pi*fi0*t+1) +3*sin(3*2*pi*fi0*t+3)+0.8*sin(5*2*pi*fi0*t+2)+0.02*sin(7*2*pi*fi0*t+1) +2*rand(size(t)); % the output is a sum of odd harmonics
x=[x1 x2];
%============ Reconstruction begins ============%
fi=sinfapm(x(:,1),fs); % the frequency is obtained by sinewave fitting.
Nr_cycles=floor(fi*length(x)/fs); % number of comlete cycles
N=round(Nr_cycles*fs/fi); % cutoff length for coherent sampling
x=x(1:N,:);
t=t(1:N); % cut the time axis too
Xraw=fft(x(:,2));
phase=angle(Xraw);
freq_bin=0:N-1; % frequency bin
freq_res=fs/N; % frequency resolution
freq_i=freq_bin*freq_res; % frequency axis
cutoff=ceil(N/2); % cut the first half of FFT
Xraw=2*Xraw(1:cutoff)/N;
phase=phase(1:cutoff);
max_nr_harm=7; % in practice the max number of harmonics is determined by the Nyquist frequency
I_n=zeros((max_nr_harm+1)/2,1); % amplitudes of the odd harmonics
delta_n=zeros((max_nr_harm+1)/2,1); % phase angles of the odd harmonics
% extract the harmonic information
for n=1:(max_nr_harm+1)/2
I_n(n)=abs(Xraw(round((2*n-1)*fi/freq_res)+1));
delta_n(n)=phase(round((2*n-1)*fi/freq_res)+1)+0.5*pi; % for the phase angles to be correct 0.5*pi must be added. The reason is unknown.
end
DC=abs(Xraw(1))*cos(angle(Xraw(1)))/2; % Normalization of the DC component need not multiply with 2, so the division by 2. The cosine factor determine the sign of DC offset.
fs2=300; % reconstruct the signal at a higher sampling frequency
t2=0:1/fs2:t(end);
t2=t2';
reconst=zeros(length(t2),2);
reconst(:,1)=sin(2*pi*fi*t2); % the input signal at high sampling frequency
% reconstruct the output signal by summing up the sine components
for n=1:2:max_nr_harm
reconst(:,2)=reconst(:,2)+I_n((n+1)/2)*sin(n*2*pi*fi*t2+delta_n((n+1)/2));
end
reconst(:,2)=reconst(:,2)+DC; % DC shift
plot(x(:,1),x(:,2),'k',reconst(:,1),reconst(:,2),'r')
为了展示，该代码生成一个示范信号。这个示范信号的频率、取样频率和长度都尽可能的任意（频率是个除不尽的小数，取样频率是个质数，长度是12345），也增加了比较明显的噪音。波形重构时，频率是通过对输入信号进行正弦波拟合来获得的。这一步对原信号的噪音还是十分敏感的。如果在以上代码给原信号加上高于千分之一的噪音，整个代码的结果就完全不及格了。这就要由于是后面决定谐波所处的frequency bin那步对信号频率的准确性十分敏感。
有了信号的频率，就可以进行相干取样、FFT……一切都按正常进行。
以上代码本身会给出一个结果比较。我用以上方法对一个真实的LAOS信号进行重构，得到的结果如下(黑的是实验结果，红的是重构曲线）：
实际信号示范
我在第7届青年流变学家研讨会做的一些记录，陆续做一些总结。由于时间有限我写我自己才看得懂的东西，不多做解释。
会议的报告表现出几个比较集中的话题。有些话题我不太熟悉，幸亏导师一直坐我旁边。任何不懂的东西都可以马上问他。有时导师会叫我直接向演讲人提问。
线性粘弹性理论仍然是重要的流变学手段，特别是对于高分子体系，包括熔体和溶液。这些体系通过时温叠加（TTS），可获得“完整”的松弛时间谱，然后再通过松弛时间谱中显示出来的快慢模式来获取分子运动的信息。从G’、G” vs ω曲线获得松弛时间分布的过程是需要非线性拟合。AR-G2带的Rheology Advantage软件就有此功能。但是大家还是一直提出各种算法来做这件事情，就是因为在数学上从动态模量获得连续松弛时间谱H(τ)需要解线性粘弹性的那个Voltera积分方程。解这个方程是不适定（ill-posed）的反问题。正则化（regularization）是一个比较常规的解决方法。基于正则化的方法有很多报道。K. Cho等提出新的算法不需要使用regularization。
在使用K-BKZ本构方程时，记忆函数往往通过线性粘弹性的实验来确定。通过TTS得到了G’、G”主曲曲线之后，是要用一个有限模式数的Maxwell方程去拟合，得到一个离散的松弛时间谱，才能把记乙函数G(t)写成一个求和式。之所以要用离散谱，是因为各种高分子体系的连续时间谱往往不是一个特定的函数，只能近似写成一个truncated的级数。K. Cho又提出了能够自动判断需要几个模式的Maxwell方程来拟合的算法，不用人为规定。这个算法是基于小波分析的。我去年有段时间也对小波分析方法感兴趣，主要是为了看能不能应用在时间依赖性样品的LAOS测试中。
凝胶化过程有很多种机理。Flory的那种是双官能团和多官能团单体缩取反应形成的，完全不考虑键的形成机率（适用于共价键）的话，可能用percolation模型（因为这个模型只描述连通性，推广到kinetics的话需要增加假设）; 如果交联键是物理键，那就要考虑键寿命，有点类似触变性流体模型中的结构理论。无论如何，以上这些凝胶化过程，都可以用聚集模型来描述，会形成分形聚集体，网络有在一定空间尺度范围内具有自相似性。还有一种凝胶化，是向高分子溶液添加交联剂形成的，例如向聚丙烯酸溶液中添加三价铁离子。由于这种凝胶并不是“生长”出来的，而是在把现成的线条粘起来，因此不一定是从聚集体发展起来的。假如是交联前是很稀的高分子溶液，高分子线团之间没有重叠和缠结，那么添加了交联剂之后就有可能是使线团与线团之间发生聚集。这时线团有点类似胶体体系中的悬浮粒子，但却不能当成硬球，甚至不能当成普通软球，因为高分子线团之间有一定的重叠能力，算是“超软球”。例如，交联发生在线两线团之间，两线团变成一个线团，自然就要按照新的线团的尺寸重新取一个最佳构象，当然，交联反应快而多的话，可能构象会来不及调整，真的就类似胶体聚集体的增长，而形成一定尺度的分形；交联反应足够慢的话，线团与线团之间是coalescence，因此没有分形结构。假如交联前是浓溶液，高分子已经形成了缠结网络了，那么交联反应的发展跟高分子的扩散速度相比当然是快好几个数量级的。因此这种凝胶其实不过就是对原缠结网络的一种固定。原缠结网络不是分形结构，交联之后自然也不会是分形结构。
还有很多各种各样的水凝胶。例如，两亲性嵌段共聚物胶束之间发生桥联、PNIPAM水溶液相变、虫状胶束、多肽自组装、超分子水凝胶等等。不同凝胶化机理，kinetics是不一样的，形成的网络结构差别也很大，要用不同的模型来解释。复旦大学Sanling Lu同学报告的凝胶体系是chitosan/β-glycerol phosphate水凝胶。向chitosan水溶液添加β-glycerol phosphate之后，就有温敏性的凝胶化行为。估计是后者通过疏水作用把chitosan交联了？从她给出的显微镜照片来看，这种水凝胶形成的网络跟colloid凝胶很类似，是基于聚集体的网络，表现为结构的空间分布很不均匀。因此其流变学性质应该是类似colloid体系的。大幅振荡剪切下应该是屈服流动，以耗散为主，所以不是什么“应变硬化”。Ewoldt和McKinley那篇文章误导了不少人。他们使用的样品，在LAOS下弹性为主，他们定义的“应变硬化”指标，不是什么样品都适用的。但是至少我就看到好几个做LAOS的人把“应变硬化”指标硬套到完全不可能的体系上去，得出不符合实际的结果。
长支化（LCB）聚合物的研究经常讨论热流变复杂性。聚合物熔体的时温叠加（TTS）的水平平移因子，可以用来计算Arrhenius的活化能Ea。对于TTS失败的体系（即热流变复杂体系），从这次会议知道原来也有很多人提出各种奇技淫巧来计算活化能。按照最一般的知识，热流变复杂性是由于样品中的松弛模式的温度依赖性不同导致的。温度依赖性就是Ea。因此，正确处理热流变复杂性样品数据的方法，应该是拿Ea对松弛模式作图。我看到有些文献的做法，是单独用G’、G”、δ、η0等曲线来组成主曲线，获得不同的活化能。这些活化能的意义还是很模糊。可是，要对松弛模式作图，首先要获得松弛时间谱。对于一个叠不出主曲线来的样品，如何从各个温度的结果得出松弛时间谱也是个难题。比较直观的做法就是，整个主曲线叠不好，我就分成至少能叠得好的几截。在频率坐标上分的这几截，每截应该对应着一个平均松弛时间的。因此每截上的主曲线对应的Ea就是这个松弛时间下的活化能，作出图来就能看到活化能对松弛时间的依赖关系。
以上讨论仅限于Arrhenius关系适用的情况。还有一种热流变复杂性是由于接近玻璃态导致的，也就是主曲线在高频部分叠得不好。这种情况Plazek和Ngai等人进行了很多研究。