前天，一個來自河北的研究生小薇（化名），提出了一個令人頭疼的繪圖問題：根據(jù)充放電曲線，怎樣用Origin的“微分”數(shù)學(xué)方法繪制“微分容量曲線”。

譚編之前讀碩士是做鋰電池材料的，不記得是藍電、輸力強還是Autolab哪一臺電化學(xué)工作站好像可以直接測試微分容量曲線。

在充放電曲線中V-t曲線包含豐富的電極過程信息，但這些信息一般很小，不容易表現(xiàn)出來，例如小薇提問的充放電曲線（如下圖）幾乎看不出充放電平臺的電位。

如果利用儀器測量是可以直接獲得微分容量曲線的，只是要注意：

1. 測試恒電流一定要盡量小，可以忽略極化的影響，一般采用0.04C倍率（即1/25 C）。在電流足夠低的情況下，可以獲得比循環(huán)伏安曲線更加豐富的反應(yīng)信息。

2. 數(shù)據(jù)采集要求：選擇電壓間隔采點，一般ΔV=10~50 mV。最好使用高精度的安捷倫、歐姆龍等儀器采集電壓（可精確到小數(shù)點后6位）。充電柜、藍電等儀器采集電壓數(shù)據(jù)誤差較大，作出的曲線毛刺多（噪音大）。

藍電儀器可以對測試的充放電曲線進行處理，得到微分容量曲線，其基本原理可以用下圖解釋。

經(jīng)過處理的微分容量曲線如下圖所示，是不是跟循環(huán)伏安曲線非常像，但比循環(huán)伏安曲線更多氧化還原峰，顯示的信息更豐富。

依據(jù)上述原理，怎樣從Originlab軟件對充放電曲線進行“微分（Differentiate）”處理（如下圖），在彈窗中直接點擊確定，一般情況下即可完成對某曲線的微分。

但是，實際上，因為充放電曲線的測試精度較低，曲線并非數(shù)學(xué)函數(shù)平滑的曲線，往往在利用Origin進行微分處理時，會出錯。

例如，假設(shè)相鄰2個數(shù)據(jù)點的Δy和Δx，有可能因為數(shù)據(jù)相等導(dǎo)致Δy或Δx→0甚至等于0，這就導(dǎo)致某一點的微分結(jié)果趨近于0或等于0（如下圖瀑布狀豎線就是逼近0），或者出現(xiàn)“被0除”的邏輯錯誤（如下圖表格中紅線框無數(shù)據(jù)，導(dǎo)致繪圖中上端整齊的“天花板”，這里200以上無數(shù)據(jù)，都是因為邏輯錯誤而“缺失”了數(shù)據(jù)）。

其實，上圖的微分結(jié)果已經(jīng)“初具雛形”（綠色畫線的形狀），譚編嘗試了很多次（就像小薇的苦惱困惑一樣）都是這樣一個結(jié)果，為什么會醬紫？

譚編初步分析的結(jié)果，還是因為充放電曲線的數(shù)據(jù)精度問題導(dǎo)致的。我請教了廣財大的我哥（數(shù)學(xué)教師），說對曲線求導(dǎo)時求導(dǎo)步長越小越精確，一般情況，相隔4個點的數(shù)據(jù)求導(dǎo)就能滿足了。

于是，譚編畫了一個示意圖如下，很理想的狀態(tài)了，相鄰數(shù)據(jù)點之間是平滑的。假設(shè)每隔2個點提取一個需要操作的數(shù)據(jù)點（即研究相鄰4點之間的偏差，它們之間的2個點不考慮），得到(x1,y1)、(x2,y2)、(x3,y3)、(x4,y4)四個操作點。

這幾個操作點兩兩之間求差：ΔV=x2-x1，ΔQ=y2-y1。

根據(jù)曲線求導(dǎo)的結(jié)果是曲線上某點的斜率（曲線上該點的切線斜率，表示該點在y方向上增量變化的加速度），則該點的斜率：

y=sinθ=ΔQ/ΔV

如果在曲線上這種操作點的選取步長較小，那么操作點之間的連線越接近真實曲線（如上圖）上點（x,y）的切線斜率，這個點（x,y）位于這兩個相鄰操作點的中點，因此，該點（x,y）的x取操作點橫坐標(biāo)的平均值：

x=(x1+x2)/2

因此，該方法暫且命名為“平均求導(dǎo)法”。曲線上某點的微分求導(dǎo)公式如下

x=(x1+x2)/2

y=sinθ=ΔQ/ΔV

根據(jù)上述公式，我們將充放電曲線的數(shù)據(jù)拷貝到Excel中進行處理。由于實際測試的實驗數(shù)據(jù)噪聲波動大，所以取間隔數(shù)d較小時，會出現(xiàn)因為相鄰數(shù)據(jù)點的大幅波動導(dǎo)致前述分析的情況。

譚編嘗試了d=10、20、30、40、100、200的間隔數(shù)（微分步長）在Excel中自動運算生成兩列新的求導(dǎo)數(shù)據(jù)。具體以d=10為例，如下圖所示，在D列第11行（其實是第10個數(shù)據(jù)）輸入“=(B11+B2)/2”，回車即可創(chuàng)建x值，點擊新建的這個x值，然后雙擊該x值單元格的右下角“+”，實現(xiàn)自動向下填充所有數(shù)據(jù)。

這里，Excel表將自動將其下方單元格的計算對象的行號加1，例如D12單元格，其計算公式為“=(B12+B3)/2”，以此類推。

同樣的方法，在求dQ/dV時，第E11單元格的計算公式為“=(C11-C2)/(B11-B2)”，雙擊向下填充即可。

按照相同的方法，計算了d的其他取值情況下的微分求導(dǎo)數(shù)據(jù)，并對這些結(jié)果繪圖如下。虛線為原始充電曲線，不同顏色的實線為不同微分步長d條件下得到的微分容量曲線。

為了比較不同d取值對微分容量峰位置的影響，將上圖利用origin繪制Stack疊圖，比較研究發(fā)現(xiàn)，當(dāng)d取值越大，得到的微分容量曲線的毛刺越來越小，當(dāng)d=200時，得到的微分容量曲線近乎完美。按道理，微分步長越大越不準(zhǔn)確，但是d=10和d=200的峰位置僅略有偏差，d越大，峰值向低壓方向漂移，漂移量ΔE=-0.019 3 V，僅減小了3.23%。因此，這種偏差可以容忍的。

上述神操作，純屬譚編瞎想獨創(chuàng)的方法，連自己都質(zhì)疑這種方法的科學(xué)性，這兩天一直不敢推出來，怕誤導(dǎo)大家，但是確實解決了頭痛的問題。

例如下圖的思路，是譚編最開始解決問題的一個思路（其實Excel自動填充公式時，就是這種效果。哈哈?。褪窍旅鎺讉€相鄰數(shù)據(jù)點（不是操作點哦），用第4點減去第1點求得第4點的斜率，用第5點減去第2點求得第5點的斜率，以此類推。這種方法暫且命名為“交叉求導(dǎo)法”。

再繪制斜率k-V曲線就能得到微分容量曲線，譚編對比研究“平均求導(dǎo)法”和“交叉求導(dǎo)法”在d=200條件下，得到的微分容量曲線的結(jié)果（如下圖）。結(jié)果非常相近，這次漂移是“正”的，漂移量ΔE=0.018 1 V，剛好跟前面d=100相對于d=10的“負(fù)”漂移相抵消。從某種程度上說，譚編獨創(chuàng)的“交叉求導(dǎo)法”的精度也能滿足我們對于“狂噪”實驗數(shù)據(jù)曲線的微分求導(dǎo)。

仔細(xì)斟酌研究這兩種方法發(fā)現(xiàn)，其實，微分步長是1，大家有沒有發(fā)現(xiàn)呢？只不過譚編的這兩種方法是跳躍采點，但每一個點都是連續(xù)求斜率的，因此，譚編的這兩種求導(dǎo)方法的微分步長就是1。