JP7026972B2 - インピーダンススペクトルデータを用いた解析処理方法、インピーダンススペクトルデータ解析処理システム、および、インピーダンススペクトル解析処理プログラム - Google Patents

インピーダンススペクトルデータを用いた解析処理方法、インピーダンススペクトルデータ解析処理システム、および、インピーダンススペクトル解析処理プログラム Download PDF

Info

Publication number
JP7026972B2
JP7026972B2 JP2020532307A JP2020532307A JP7026972B2 JP 7026972 B2 JP7026972 B2 JP 7026972B2 JP 2020532307 A JP2020532307 A JP 2020532307A JP 2020532307 A JP2020532307 A JP 2020532307A JP 7026972 B2 JP7026972 B2 JP 7026972B2
Authority
JP
Japan
Prior art keywords
logarithmic
relaxation time
measured
impedance
reset
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
JP2020532307A
Other languages
English (en)
Other versions
JPWO2020022124A1 (ja
Inventor
清 小林
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
National Institute for Materials Science
Original Assignee
National Institute for Materials Science
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by National Institute for Materials Science filed Critical National Institute for Materials Science
Publication of JPWO2020022124A1 publication Critical patent/JPWO2020022124A1/ja
Application granted granted Critical
Publication of JP7026972B2 publication Critical patent/JP7026972B2/ja
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01R—MEASURING ELECTRIC VARIABLES; MEASURING MAGNETIC VARIABLES
    • G01R31/00—Arrangements for testing electric properties; Arrangements for locating electric faults; Arrangements for electrical testing characterised by what is being tested not provided for elsewhere
    • G01R31/36—Arrangements for testing, measuring or monitoring the electrical condition of accumulators or electric batteries, e.g. capacity or state of charge [SoC]
    • G01R31/389—Measuring internal impedance, internal conductance or related variables
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01R—MEASURING ELECTRIC VARIABLES; MEASURING MAGNETIC VARIABLES
    • G01R27/00—Arrangements for measuring resistance, reactance, impedance, or electric characteristics derived therefrom
    • G01R27/28—Measuring attenuation, gain, phase shift or derived characteristics of electric four pole networks, i.e. two-port networks; Measuring transient response
    • H—ELECTRICITY
    • H01—ELECTRIC ELEMENTS
    • H01M—PROCESSES OR MEANS, e.g. BATTERIES, FOR THE DIRECT CONVERSION OF CHEMICAL ENERGY INTO ELECTRICAL ENERGY
    • H01M10/00—Secondary cells; Manufacture thereof
    • H01M10/04—Construction or manufacture in general
    • H01M10/0422—Cells or battery with cylindrical casing
    • H—ELECTRICITY
    • H01—ELECTRIC ELEMENTS
    • H01M—PROCESSES OR MEANS, e.g. BATTERIES, FOR THE DIRECT CONVERSION OF CHEMICAL ENERGY INTO ELECTRICAL ENERGY
    • H01M10/00—Secondary cells; Manufacture thereof
    • H01M10/05—Accumulators with non-aqueous electrolyte
    • H01M10/052—Li-accumulators
    • H—ELECTRICITY
    • H01—ELECTRIC ELEMENTS
    • H01M—PROCESSES OR MEANS, e.g. BATTERIES, FOR THE DIRECT CONVERSION OF CHEMICAL ENERGY INTO ELECTRICAL ENERGY
    • H01M10/00—Secondary cells; Manufacture thereof
    • H01M10/42—Methods or arrangements for servicing or maintenance of secondary cells or secondary half-cells
    • H01M10/4285—Testing apparatus
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01R—MEASURING ELECTRIC VARIABLES; MEASURING MAGNETIC VARIABLES
    • G01R31/00—Arrangements for testing electric properties; Arrangements for locating electric faults; Arrangements for electrical testing characterised by what is being tested not provided for elsewhere
    • G01R31/36—Arrangements for testing, measuring or monitoring the electrical condition of accumulators or electric batteries, e.g. capacity or state of charge [SoC]
    • G01R31/367—Software therefor, e.g. for battery testing using modelling or look-up tables
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01R—MEASURING ELECTRIC VARIABLES; MEASURING MAGNETIC VARIABLES
    • G01R31/00—Arrangements for testing electric properties; Arrangements for locating electric faults; Arrangements for electrical testing characterised by what is being tested not provided for elsewhere
    • G01R31/50—Testing of electric apparatus, lines, cables or components for short-circuits, continuity, leakage current or incorrect line connections
    • G01R31/64—Testing of capacitors
    • H—ELECTRICITY
    • H01—ELECTRIC ELEMENTS
    • H01G—CAPACITORS; CAPACITORS, RECTIFIERS, DETECTORS, SWITCHING DEVICES, LIGHT-SENSITIVE OR TEMPERATURE-SENSITIVE DEVICES OF THE ELECTROLYTIC TYPE
    • H01G11/00—Hybrid capacitors, i.e. capacitors having different positive and negative electrodes; Electric double-layer [EDL] capacitors; Processes for the manufacture thereof or of parts thereof
    • H01G11/14—Arrangements or processes for adjusting or protecting hybrid or EDL capacitors
    • H—ELECTRICITY
    • H01—ELECTRIC ELEMENTS
    • H01G—CAPACITORS; CAPACITORS, RECTIFIERS, DETECTORS, SWITCHING DEVICES, LIGHT-SENSITIVE OR TEMPERATURE-SENSITIVE DEVICES OF THE ELECTROLYTIC TYPE
    • H01G13/00—Apparatus specially adapted for manufacturing capacitors; Processes specially adapted for manufacturing capacitors not provided for in groups H01G4/00 - H01G11/00
    • H—ELECTRICITY
    • H01—ELECTRIC ELEMENTS
    • H01M—PROCESSES OR MEANS, e.g. BATTERIES, FOR THE DIRECT CONVERSION OF CHEMICAL ENERGY INTO ELECTRICAL ENERGY
    • H01M10/00—Secondary cells; Manufacture thereof
    • H01M10/42—Methods or arrangements for servicing or maintenance of secondary cells or secondary half-cells
    • H01M10/48—Accumulators combined with arrangements for measuring, testing or indicating the condition of cells, e.g. the level or density of the electrolyte
    • H—ELECTRICITY
    • H01—ELECTRIC ELEMENTS
    • H01M—PROCESSES OR MEANS, e.g. BATTERIES, FOR THE DIRECT CONVERSION OF CHEMICAL ENERGY INTO ELECTRICAL ENERGY
    • H01M8/00—Fuel cells; Manufacture thereof
    • H01M8/04—Auxiliary arrangements, e.g. for control of pressure or for circulation of fluids
    • H01M8/04298—Processes for controlling fuel cells or fuel cell systems
    • H01M8/04313—Processes for controlling fuel cells or fuel cell systems characterised by the detection or assessment of variables; characterised by the detection or assessment of failure or abnormal function
    • H01M8/04537—Electric variables
    • Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
    • Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
    • Y02E—REDUCTION OF GREENHOUSE GAS [GHG] EMISSIONS, RELATED TO ENERGY GENERATION, TRANSMISSION OR DISTRIBUTION
    • Y02E60/00—Enabling technologies; Technologies with a potential or indirect contribution to GHG emissions mitigation
    • Y02E60/30—Hydrogen technology
    • Y02E60/50—Fuel cells
    • Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
    • Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
    • Y02P—CLIMATE CHANGE MITIGATION TECHNOLOGIES IN THE PRODUCTION OR PROCESSING OF GOODS
    • Y02P70/00—Climate change mitigation technologies in the production process for final industrial or consumer products
    • Y02P70/50—Manufacturing or production processes characterised by the final manufactured product

Landscapes

  • Engineering & Computer Science (AREA)
  • Manufacturing & Machinery (AREA)
  • Chemical & Material Sciences (AREA)
  • Chemical Kinetics & Catalysis (AREA)
  • Electrochemistry (AREA)
  • General Chemical & Material Sciences (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Measurement Of Resistance Or Impedance (AREA)
  • Secondary Cells (AREA)
  • Primary Cells (AREA)
  • Investigating Or Analyzing Materials By The Use Of Electric Means (AREA)

Description

本発明は、インピーダンススペクトルデータを用いた解析処理方法、インピーダンススペクトルデータ解析処理システム、および、インピーダンススペクトル解析処理プログラムに関する。
測定されたインピーダンススペクトルデータの解析法として、連続緩和時間分布を求める方法が知られている（例えば、非特許文献１を参照）。非特許文献１によれば、対数で等間隔に取られた周波数（ｆ）と複素インピーダンス（Ｚ）とが測定データとして得られた場合、複素インピーダンスは次式（ａ）および（ｂ）で与えられるとして解析される。式（ａ）のＲ∞とＲｋは未定値であり測定データを再現するために必要なパラメータである。
Figure 0007026972000001
ここで、ｋは対数緩和時間の配列番号、ｐは測定されたデータ点の配列番号、ｆｋおよびｆｐはｋ番目およびｐ番目の周波数、ｊは単位虚数、Ｎは測定データ点数、Ｒ∞は高周波限抵抗、τｋはｋ番目の緩和時間、Ｒｋはτｋにおける抵抗値、｜Δｌｎτ｜は自然対数緩和時間の間隔の絶対値である。
しかしながら、非特許文献１によれば、一次電池、二次電池、キャパシタ、腐食系のインピーダンスのような周波数の低下によってインピーダンス絶対値が一定の値に収束しない、または、インピーダンス絶対値が一定の値に収束していないデータに対して、式（ａ）および式（ｂ）を適用し得られるＲｋ値セットは、測定されたインピーダンススペクトルデータを再現できない場合がある。さらに、上述の試料はインダクタ成分を含むため、インピーダンス虚数項が正値を取る場合にも式（ａ）および式（ｂ）は適用できない。
Ｐ．Ｂｒａｕｎら，Ｊ． Ｅｌｅｃｔｒｏｃｅｒａｍ．，３８，１５７－１６７（２０１７）
以上から、本発明の課題は、測定されたインピーダンススペクトルデータを用いた解析処理方法、インピーダンススペクトルデータ解析処理システム、および、インピーダンススペクトルデータ解析処理プログラムを提供することである。本発明のさらなる課題は、測定されたインピーダンススペクトルデータに対して緩和時間分布解析を可能にする解析処理方法、インピーダンススペクトルデータ解析処理システム、および、インピーダンススペクトルデータ解析処理プログラムを提供することである。
本発明による測定されたインピーダンススペクトルデータを用いて解析処理する方法は、前記測定されたインピーダンススペクトルデータに基づいて、測定された対数周波数に対応する対数緩和時間の最大値および／または最小値を設定し直すか否かを判定するステップと、前記判定するステップで設定し直すと判定した場合、前記対数緩和時間の最大値および／または最小値を設定し直すステップと、前記判定するステップで設定し直さないと判定した場合には前記測定された対数周波数に対応する対数緩和時間の範囲、または、前記判定するステップで設定し直すと判定した場合には前記設定し直された対数緩和時間の範囲で等間隔に分割する分割数を設定するステップと、正則化最小二乗法を適用し、次式（Ａ）
Figure 0007026972000002
（ｐは測定されたデータ点の配列番号、ｌは対数緩和時間の配列番号、ｆｐはｐ番目の周波数、ｊは単位虚数、Ｍは対数緩和時間の分割数、Ｔはキャパシタンス（Ｃ）の逆数、Ｌはインダクタンス、Ｒ∞は高周波限抵抗、τｌはｌ番目の緩和時間、Ｒｌはτｌにおける抵抗値、｜Δｌｎτ｜は自然対数緩和時間の間隔の絶対値）
を満たすように、前記対数緩和時間の範囲および設定した前記分割数について、パラメータＲ∞、Ｔ、ＬおよびＲｌを解析するステップとを包含し、これにより上記課題を解決する。
前記判定するステップは、前記測定されたインピーダンススペクトルデータが、（１）インピーダンス実数項と対数周波数との関係において、インピーダンス実数項が、測定された対数周波数の範囲のうち低周波数側および高周波数側でプラトーを有さないか、（２）インピーダンス虚数項と対数周波数との関係において、インピーダンス虚数項が、測定された対数周波数の範囲のうち低周波数側および高周波数側でピークを有さないか、または、（３）インピーダンスプロットにおいて円弧を有さない場合に、前記解析するステップに用いる対数緩和時間の最大値および／または最小値を設定し直すと判定し、前記設定し直すステップは、（１’）前記インピーダンス実数項と対数周波数との関係において、前記インピーダンス実数項が、前記測定された対数周波数の範囲を超えて低周波数側および高周波数側でプラトーとなるように前記解析するステップに用いる対数緩和時間の最大値および／または最小値を設定し直すか、（２’）前記インピーダンス虚数項と対数周波数との関係において、前記インピーダンス虚数項が、前記測定された対数周波数の範囲を超えて低周波数側および高周波数側でピークを有するように前記解析するステップに用いる対数緩和時間の最大値および／または最小値を設定し直すか、（３’）前記インピーダンスプロットにおいて円弧を有するように前記解析するステップに用いる対数緩和時間の最大値および／または最小値を設定し直してもよい。
前記判定するステップは、前記測定された対数周波数の範囲について、次式（ａ）および次式（ｂ）
Figure 0007026972000003
（ｋは対数緩和時間の配列番号、ｐは測定されたデータ点の配列番号、ｆｋおよびｆｐはｋ番目およびｐ番目の周波数、ｊは単位虚数、Ｎは測定データ点数、Ｒ∞は高周波限抵抗、τｋはｋ番目の緩和時間、Ｒｋはτｋにおける抵抗値、｜Δｌｎτ｜は自然対数緩和時間の間隔の絶対値）
を用いて、自然対数緩和時間分布抵抗Ｒｋの対数緩和時間依存性を求め、前記測定された対数周波数に相当する対数緩和時間の範囲に前記Ｒｋがピークを有さない場合に、前記解析するステップに用いる対数緩和時間の最大値および／または最小値を設定し直すと判定し、前記設定し直すステップは、前記測定された対数周波数に相当する対数緩和時間の範囲を超えて前記Ｒｋがピークを有するように、前記解析するステップに用いる対数緩和時間の最大値および／または最小値を設定し直してもよい。
前記解析するステップにおいて、前記式（Ａ）における前記キャパシタのインピーダンスは、次式（Ｂ）
Ｚ（ｆ）＝Ｔ／（ｊ２πｆ） ・・・（Ｂ）
で与えられてもよい。
前記解析するステップは、正則化非線形最適化法を用いて算出してもよい。
前記正則化非線形最適化法は、Ｔｉｋｈｏｎｏｖ正則化法、Ｌａｓｓｏ正則化法およびＥｌａｓｔｉｃ Ｎｅｔ正則化法からなる群から選択される正則化法を採用してもよい。
あるいは、前記正則化非線形最適化法は、Ｌｅｖｅｎｂｅｒｇ－Ｍａｒｑｕａｒｄｔ法、シンプレックス法、および、共役勾配法からなる群から選択される非線形最適化法を採用してもよい。
前記解析するステップで得られた前記パラメータを、前記式（Ａ）に代入し、理論インピーダンススペクトルを求めるステップと、前記求めた理論インピーダンススペクトルと前記測定されたインピーダンススペクトルとを表示し、比較し、一致するか否かを判定するステップとをさらに包含してもよい。
前記表示し、比較し、一致するか否かを判定するステップが、一致しないと判定した場合、前記設定し直すステップ以降を繰り返してもよい。
前記測定されたインピーダンススペクトルデータは、一次電池、二次電池、燃料電池、キャパシタ材料、防蝕処理を施した金属およびセラミックス材料からなる群から選択される試料のインピーダンススペクトルデータであってもよい。
また、本発明による測定されたインピーダンススペクトルデータを用いたインピーダンススペクトルデータ解析処理システムは、制御装置と、表示装置とを備え、前記制御装置は、前記測定されたインピーダンススペクトルデータに基づいて、測定された対数周波数に対応する対数緩和時間の最大値および／または最小値を設定し直すか否かを判定する対数緩和時間判定部と、前記対数緩和時間判定部が設定し直すと判定した場合、前記測定された対数周波数に対応する対数緩和時間の最大値および／または最小値を設定し直す対数緩和時間設定部と、前記対数緩和時間判定部が設定し直さないと判定した場合には前記測定された対数周波数に対応する対数緩和時間の範囲、または、前記対数緩和時間判定部が設定し直すと判定し場合には前記対数緩和時間設定部が設定し直した対数緩和時間の範囲で等間隔に分割する分割数を設定する分割数設定部と、正則化最小二乗法を適用し、次式（Ａ）
Figure 0007026972000004
（ｐは測定されたデータ点の配列番号、ｌは対数緩和時間の配列番号、ｆｐはｐ番目の周波数、ｊは単位虚数、Ｍは対数緩和時間の分割数、Ｔはキャパシタンス（Ｃ）の逆数、Ｌはインダクタンス、Ｒ∞は高周波限抵抗、τｌはｌ番目の緩和時間、Ｒｌはτｌにおける抵抗値、｜Δｌｎτ｜は自然対数緩和時間の間隔の絶対値）
を満たすように、前記対数緩和時間の範囲および設定した前記分割数について、パラメータＲ∞、Ｔ、ＬおよびＲｌを解析するパラメータ解析部とをさらに備え、前記表示装置は、前記制御装置における解析結果を表示し、これにより上記課題を解決する。
前記制御装置の前記対数緩和時間判定部は、前記測定されたインピーダンススペクトルデータが、（１）インピーダンス実数項と対数周波数との関係において、インピーダンス実数項が、測定された対数周波数の範囲のうち低周波数側および高周波数側でプラトーを有さないか、（２）インピーダンス虚数項と対数周波数との関係において、インピーダンス虚数項が、測定された対数周波数の範囲のうち低周波数側および高周波数側でピークを有さないか、または、（３）インピーダンスプロットにおいて円弧を有さない場合に、前記パラメータ解析部で用いる対数緩和時間の最大値および／または最小値を設定し直すと判定し、前記制御装置の前記対数緩和時間設定部は、（１’）前記インピーダンス実数項と対数周波数との関係において、前記インピーダンス実数項が、前記測定された対数周波数の範囲を超えて低周波数側および高周波数側でプラトーとなるように、前記パラメータ解析部で用いる対数緩和時間の最大値および／または最小値を設定し直すか、（２’）前記インピーダンス虚数項と対数周波数との関係において、前記インピーダンス虚数項が、前記測定された対数周波数の範囲を超えて低周波数側および高周波数側でピークを有するように、前記パラメータ解析部で用いる対数緩和時間の最大値および／または最小値を設定し直すか、（３’）前記インピーダンスプロットにおいて円弧を有するように前記パラメータ解析部で用いる対数緩和時間の最大値および／または最小値を設定し直してもよい。
前記制御装置の前記対数緩和時間判定部は、前記測定された対数周波数の範囲について、次式（ａ）および次式（ｂ）
Figure 0007026972000005
（ｋは対数緩和時間の配列番号、ｐは測定されたデータ点の配列番号、ｆｋおよびｆｐはｋ番目およびｐ番目の周波数、ｊは単位虚数、Ｎは測定データ点数、Ｒ∞は高周波限抵抗、τｋはｋ番目の緩和時間、Ｒｋはτｋにおける抵抗値、｜Δｌｎτ｜は自然対数緩和時間の間隔の絶対値）を用いて、自然対数緩和時間分布抵抗Ｒｋの対数緩和時間依存性を求め、前記測定された対数周波数の範囲に前記Ｒｋがピークを有さない場合に、前記パラメータ解析部で用いる対数緩和時間の最大値および／または最小値を設定し直すと判定し、前記制御装置の前記対数緩和時間設定部は、前記測定された対数周波数に相当する対数緩和時間の範囲を超えて前記Ｒｋがピークを有するように前記パラメータ解析部で用いる対数緩和時間の最大値および／または最小値を設定し直してもよい。
前記制御装置は、オペレータが、前記表示装置に表示された前記測定されたインピーダンススペクトルデータから目視にて前記（１）～（３）を判定し、前記（１’）～（３’）を手動にて設定し直す入力端末を備えてもよい。
前記制御装置は、オペレータが、前記表示装置に表示された前記自然対数緩和時間分布抵抗Ｒｋの対数緩和時間依存性から目視にて前記測定された対数周波数の範囲に前記Ｒｋがピークを有するか否かを判定し、前記測定された対数周波数の範囲を超えて前記Ｒｋがピークを有するように前記パラメータ解析部で用いる対数緩和時間の最大値および／または最小値を手動にて設定し直す入力端末を備えてもよい。
前記式（Ａ）において、前記キャパシタのインピーダンスは、次式（Ｂ）
Ｚ（ｆ）＝Ｔ／（ｊ２πｆ） ・・・（Ｂ）
で与えられてもよい。
前記制御装置は、得られた前記パラメータを、前記式（Ａ）に代入し、理論インピーダンススペクトルを求め、前記求めた理論インピーダンススペクトルと前記測定されたインピーダンススペクトルとを前記表示装置に表示し、比較し、前記理論インピーダンススペクトルと前記測定されたインピーダンススペクトルとが一致するか否かを判定する判定部をさらに備えてもよい。
また、本発明によるインピーダンススペクトルデータ解析処理プログラムは、コンピュータに、測定されたインピーダンススペクトルデータに基づいて、測定された対数周波数に対応する対数緩和時間の最大値および／または最小値を設定し直すか否かを判定するステップと、前記判定するステップで設定し直すと判定した場合、前記対数緩和時間の最大値および／または最小値を設定し直すステップと、前記判定するステップで設定し直さないと判定した場合には前記測定された対数周波数に対応する対数緩和時間の範囲、または、前記判定するステップで設定し直すと判定した場合には前記設定し直された対数緩和時間の範囲で等間隔に分割する数を設定するステップと、正則化最小二乗法を適用し、次式（Ａ）
Figure 0007026972000006
（ｐは測定されたデータ点の配列番号、ｌは対数緩和時間の配列番号、ｆｐはｐ番目の周波数、ｊは単位虚数、Ｍは対数緩和時間の分割数、Ｔはキャパシタンス（Ｃ）の逆数、Ｌはインダクタンス、Ｒ∞は高周波限抵抗、τｌはｌ番目の緩和時間、Ｒｌはτｌにおける抵抗値、｜Δｌｎτ｜は自然対数緩和時間の間隔の絶対値）
を満たすように、前記対数緩和時間の範囲および設定した分割数について、パラメータＲ∞、Ｔ、ＬおよびＲｌを解析するステップとを実行させ、これにより上記課題を解決する。
前記判定するステップは、前記測定されたインピーダンススペクトルデータが、（１）インピーダンス実数項と対数周波数との関係において、インピーダンス実数項が、測定された対数周波数の範囲のうち低周波数側および高周波数側でプラトーを有さないか、（２）インピーダンス虚数項と対数周波数との関係において、インピーダンス虚数項が、測定された対数周波数の範囲のうち低周波数側および高周波数側でピークを有さないか、または（３）インピーダンスプロットにおいて円弧を有さない場合に、前記解析するステップに用いる対数緩和時間の最大値および／または最小値を設定し直すと判定し、前記設定し直すステップは、（１’）前記インピーダンス実数項と対数周波数との関係において、前記インピーダンス実数項が、前記測定された対数周波数の範囲を超えて低周波数側および高周波数側でプラトーとなるように前記解析するステップに用いる対数緩和時間の最大値および／または最小値を設定し直すか、（２’）前記インピーダンス虚数項と対数周波数との関係において、前記インピーダンス虚数項が、前記測定された対数周波数の範囲を超えて低周波数側および高周波数側でピークを有するように前記解析するステップに用いる対数緩和時間の最大値および／または最小値を設定し直すか、（３’）前記インピーダンスプロットにおいて円弧を有するように前記解析するステップに用いる対数緩和時間の最大値および／または最小値を設定し直してもよい。
前記判定するステップは、測定された対数周波数の範囲について、次式（ａ）および次式（ｂ）
Figure 0007026972000007
（ｋは対数緩和時間の配列番号、ｐは測定されたデータ点の配列番号、ｆｋおよびｆｐはｋ番目およびｐ番目の周波数、ｊは単位虚数、Ｎは測定データ点数、Ｒ∞は高周波限抵抗、τｋはｋ番目の緩和時間、Ｒｋはτｋにおける抵抗値、｜Δｌｎτ｜は自然対数緩和時間の間隔の絶対値）
を用いて、自然対数緩和時間分布抵抗Ｒｋの対数緩和時間依存性を求め、前記測定された対数周波数の範囲に前記Ｒｋがピークを有さない場合に、前記解析するステップに用いる対数緩和時間の最大値および／または最小値を設定し直すと判定し、前記設定し直すステップは、前記測定された対数周波数に相当する対数緩和時間の範囲を超えて前記Ｒｋがピークを有するように前記解析するステップに用いる対数緩和時間の最大値および／または最小値を設定し直してもよい。
前記解析するステップで得られた前記パラメータを、前記式（Ａ）に代入し、理論インピーダンススペクトルを求めるステップと、前記求めた理論インピーダンススペクトルと前記測定されたインピーダンススペクトルとを表示装置に表示し、比較し、一致するか否かを判定するステップとをさらに実行させてもよい。
本発明の測定されたインピーダンススペクトルデータを用いて解析処理する方法は、測定されたインピーダンススペクトルデータに基づいて、計算に用いる対数緩和時間の最大値および／または最小値を設定し直すか否かを判定するステップと、設定した対数緩和時間の範囲で等間隔に分割する数を設定するステップと、正則化最小二乗法を適用し、上述の式（Ａ）を満たすように、設定した対数緩和時間の範囲について、パラメータＲ∞、Ｔ、ＬおよびＲｌを解析するステップとを包含する。
本発明者は、キャパシタ成分およびインダクタ成分を考慮した新たな式（Ａ）を採用することにより、一次電池、二次電池、腐食系、キャパシタのインピーダンスのように周波数の低下によってインピーダンス絶対値が一定の値に収束しない、または、インピーダンス絶対値が一定の値に収束していないデータであっても、測定されたインピーダンススペクトルデータを再現する各種パラメータを求めることができることを見出した。このようにして求めたパラメータを利用することにより、上述のインピーダンススペクトルデータに対しても緩和時間分布解析が可能となる。また、計算に用いる対数緩和時間の最大値および／または最小値を設定し直すか否かを判定し、設定した対数緩和時間領域での等間隔分割数を設定するので、上記式（Ａ）に基づいてパラメータを求めた際に、高精度に測定されたインピーダンススペクトルデータを再現できる。このような解析処理方法を採用すれば、解析処理システムおよびそれを実行するプログラムを提供できる。
図１は、本発明のインピーダンススペクトルデータを解析処理する工程を示すフローチャートである。 図２は、判定するステップを説明する模式図である。 図３は、別の判定するステップを説明する模式図である。 図４は、図１に続く本発明のインピーダンススペクトルデータを解析処理する工程を示すフローチャートである。 図５は、本発明のインピーダンススペクトルデータ解析処理システムを示すブロック図である。 図６は、比較例１で用いた電気回路を示す図である。 図７は、図６に示す電気回路を用いて１ＭＨｚから１Ｈｚの周波数で計算したインピーダンススペクトルデータの一覧である。 図８は、比較例１による図７のインピーダンススペクトルデータに基づく自然対数緩和時間分布抵抗Ｒｋの対数緩和時間依存性を示す図である。 図９は、比較例１による図８で得られた緩和時間分布解析から計算した理論インピーダンススペクトルデータと、測定されたインピーダンススペクトルデータとの比較を示す図である。 図１０は、実施例１による図７のインピーダンススペクトルデータに基づいて、対数緩和時間の範囲を再設定した自然対数緩和時間分布抵抗Ｒｌの対数緩和時間依存性を示す図である。 図１１は、実施例１による図１０で得られた緩和時間分布解析から計算した理論インピーダンススペクトルデータと、測定されたインピーダンススペクトルデータとの比較を示す図である。 図１２は、比較例２のコイン型リチウム電池のインピーダンススペクトルデータを示す図である。 図１３は、比較例２による図１２のインピーダンススペクトルデータに基づく自然対数緩和時間分布抵抗Ｒｋの対数緩和時間依存性を示す図である。 図１４は、比較例２による図１３で得られた緩和時間分布解析から計算した理論インピーダンススペクトルデータと、測定されたインピーダンススペクトルデータとの比較を示す図である。 図１５は、実施例２による図１２のインピーダンススペクトルデータに基づいて、対数緩和時間の範囲を再設定した自然対数緩和時間分布抵抗Ｒｌの対数緩和時間依存性を示す図である。 図１６は、実施例２による図１５で得られた緩和時間分布解析から計算した理論インピーダンススペクトルデータと、測定されたインピーダンススペクトルデータとの比較を示す図である。 図１７は、比較例３のＧｏｏｇｌｅ Ｉｍｐｅｄａｎｃｅ Ｓｐｅｃｔｒｏｓｃｏｐｙフォーラムにアップロードされた絶縁性酸化物のインピーダンススペクトルデータを示す図である。 図１８は、比較例３による図１７のインピーダンススペクトルデータに基づく自然対数緩和時間分布抵抗Ｒｋの対数緩和時間依存性を示す図である。 図１９は、比較例３による図１８で得られた緩和時間分布解析から計算した理論インピーダンススペクトルデータと、アップロードされたインピーダンススペクトルデータとの比較を示す図である。 図２０は、実施例３による図１７のインピーダンススペクトルデータに基づいて、対数緩和時間の範囲を再設定した自然対数緩和時間分布抵抗Ｒｌの対数緩和時間依存性を示す図である。 図２１は、実施例３による図２０で得られた緩和時間分布解析から計算した理論インピーダンススペクトルデータと、アップロードされたインピーダンススペクトルデータとの比較を示す図である。
以下、図面を参照しながら本発明の実施の形態を説明する。なお、同様の要素には同様の番号を付し、その説明を省略する。
図１は、本発明のインピーダンススペクトルデータを解析処理する工程を示すフローチャートである。
ステップＳ１１０：測定されたインピーダンススペクトルデータに基づいて、測定された対数周波数に対応する対数緩和時間の最大値および／または最小値を設定し直すか否かを判定する。後述する解析に用いる対数緩和時間の範囲の見直しを行うことにより、解析の精度を向上させることができる。なお、測定されたインピーダンススペクトルデータは、電気化学測定によって測定された外部データを用いてもよいし、その場で測定したデータであってもよい。ステップＳ１１０で設定し直すと判定した場合は、次のステップＳ１２０を行うが、設定し直さないと判定した場合はステップＳ１３０を行う。
ステップＳ１２０：測定された対数周波数に対応する対数緩和時間の最大値および／または最小値を設定し直す。
ステップＳ１３０：ステップＳ１２０で設定し直した対数緩和時間の範囲、または、ステップＳ１１０で設定し直さないと判定した場合には測定された対数周波数に対応する対数緩和時間の範囲で等間隔に分割する数を設定する。これにより、対数緩和時間の分割点数を制御できるので、解析の精度や解析速度を向上させることができる。このような分割数の設定は、オペレータ（ユーザ）が手動にて設定してもよいし、コンピュータが適宜設定するようにしてもよい。
ステップＳ１４０：正則化最小二乗法を適用し、式（Ａ）を満たすように、ステップＳ１２０で設定した対数緩和時間の範囲およびステップＳ１３０で設定した分割数について、パラメータＲ∞、Ｔ、ＬおよびＲｌを解析する。
Figure 0007026972000008
ここで、ｐは測定されたデータ点の配列番号、ｌは対数緩和時間の配列番号、ｆｐはｐ番目の周波数、ｊは単位虚数、Ｍは対数緩和時間の分割数、Ｔはキャパシタンス（Ｃ）の逆数、Ｌはインダクタンス、Ｒ∞は高周波限抵抗、τｌはｌ番目の緩和時間、Ｒｌはτｌにおける抵抗値、｜Δｌｎτ｜は自然対数緩和時間の間隔の絶対値である。
式（Ａ）は、本願発明者が新たに考案した式であり、キャパシタ成分およびインダクタ成分を考慮しているため、周波数の低下によってインピーダンス絶対値が一定の値に収束しない、または、インピーダンス絶対値が一定の値に収束していないデータに対しても、測定されたインピーダンススペクトルデータを再現する各種パラメータを求めることができる。さらに、式（Ａ）の解析に用いる対数緩和時間の範囲およびその分割数は、ステップＳ１２０およびステップＳ１３０において見直されているため、各種パラメータを高精度に求めることができる。
ここで、ステップＳ１１０における判定方法について説明する。
図２は、判定するステップを説明する模式図である。
図２の（Ａ）は、測定されたインピーダンススペクトルデータのインピーダンス実数項と対数周波数との関係を模式的に示す図である。図２の（Ｂ）は、測定されたインピーダンススペクトルデータのインピーダンス虚数項と対数周波数との関係を模式的に示す図である。図２の（Ｃ）は、測定されたインピーダンススペクトルデータのインピーダンスプロットを模式的に示す図である。図２の（Ａ）～（Ｃ）において、実線が測定されたインピーダンススペクトルデータであるとする。
図２の（Ａ）において、（１）インピーダンス実数項が、測定された対数周波数範囲において低周波数側および高周波数側でプラトー（横ばい）でない場合には、測定された対数周波数に対応する対数緩和時間の最大値および／または最小値を設定し直すと判定する。設定し直すと判定した場合には、（１’）インピーダンス実数項が測定された対数周波数の範囲を超えて低周波数側および高周波数側でプラトーとなるように、測定された対数周波数に対応する対数緩和時間の最大値および／または最小値を設定し直せばよい。
なお、本願明細書において、「プラトー（横ばい）であるか否か」とは、後述する実施例のように、インピーダンス実数項の対数周波数に対する変化を目視にて観察した際に、高周波数側および低周波数側において横ばいの領域を有するか否かを意図するが、インピーダンス実数項を対数でとった際の変化が高周波数側および低周波数側で直線的になっている（プラトーである）、直線的になっていない（プラトーでない）と判定してもよい。このような対数を用いた判定は後述する制御装置の対数緩和時間判定部が行ってもよい。
図２の（Ａ）を参照すれば、インピーダンス実数項が、低周波数側（ｆｍｉｎ）および高周波数側（ｆｍａｘ）のいずれでもプラトーを示していないため、測定された対数周波数の範囲を超えて低周波数側（ｆｍｉｎ’）および高周波数側（ｆｍａｘ’）でプラトーとなるようにし、これを対応する対数緩和時間の最大値および最小値とすればよい。なお、図２の（Ａ）では、測定された対数周波数の低周波数側および高周波数側のいずれも設定し直すと判定したが、いずれかのインピーダンス実数項がプラトーであれば一方だけ設定し直してもよい。
図２の（Ｂ）において、（２）インピーダンス虚数項が、測定された対数周波数範囲において低周波数側および高周波数側で下向きピークを有しない場合には、測定された対数周波数に対応する対数緩和時間の最大値および／または最小値を設定し直すと判定する。設定し直すと判定した場合には、（２’）インピーダンス虚数項が測定された対数周波数の範囲を超えて低周波数側および高周波数側で下向きピークを有するように、測定された対数周波数に対応する対数緩和時間の最大値および／または最小値を設定し直せばよい。
図２の（Ｂ）を参照すれば、インピーダンス虚数項が、低周波数側（ｆｍｉｎ）および高周波数側（ｆｍａｘ）のいずれでも下向きピークを有していないため、測定された対数周波数の範囲を超えて低周波数側（ｆｍｉｎ’）および高周波数側（ｆｍａｘ’）で下向きピークを有するようにし、これを対応する対数緩和時間の最大値および最小値とすればよい。なお、図２の（Ｂ）では、測定された対数周波数の低周波数側および高周波数側のいずれも設定し直すと判定したが、いずれかのインピーダンス虚数項が下向きピークを有していれば一方だけ設定し直してもよい。
図２の（Ｃ）において、（３）インピーダンスプロットにおいて円弧を有しない場合には、測定された対数周波数に対応する対数緩和時間の最大値および／または最小値を設定し直すと判定する。設定し直すと判定した場合には、（３’）インピーダンスプロットにおいて円弧を有するように、例えば、インピーダンス実数項をＺｒｅａｌ ｍｉｎからＺｒｅａｌ ｍｉｎ’へ、Ｚｒｅａｌ ｍａｘからＺｒｅａｌ ｍａｘ’へと変更し、これに対応する対数緩和時間の最大値および最小値とすればよい。なお、図２の（Ｃ）では、対数緩和時間の最大値および最小値のいずれも設定し直すと判定したが、インピーダンスプロットのいずれかが円弧になっていれば一方だけ設定し直してもよい。
図３は、別の判定するステップを説明する模式図である。
まず、測定されたインピーダンススペクトルデータを、測定された対数周波数の範囲について、次式（ａ）および次式（ｂ）を用いて、自然対数緩和時間分布抵抗Ｒｋの対数緩和時間依存性を求める。
Figure 0007026972000009
ここで、ｋは対数緩和時間の配列番号、ｐは測定されたデータ点の配列番号、ｆｋおよびｆｐはｋ番目およびｐ番目の周波数、ｊは単位虚数、Ｎは測定データ点数、Ｒ∞は高周波限抵抗、τｋはｋ番目の緩和時間、Ｒｋはτｋにおける抵抗値、｜Δｌｎτ｜は自然対数緩和時間の間隔の絶対値である。
式（ａ）および（ｂ）は、非特許文献１においても適用される公知の式である。このようにして求めた対数緩和時間分布解析において、測定された対数周波数に相当する対数緩和時間の範囲にＲｋのピークを有しない場合には、測定された対数周波数に対応する対数緩和時間の最大値および／または最小値を設定し直すと判定する。設定し直すと判定した場合には、測定された対数周波数に相当する対数緩和時間の範囲を超えてＲｋがあたかもピークを有するように、測定された対数周波数に対応する対数緩和時間の最大値および／または最小値を設定し直せばよい。
図３を参照すれば、Ｒkが対数緩和時間の両端（τｍｉｎおよびτｍａｘ）においてピークを有しないため、測定された対数周波数の範囲をこえてＲｋが対数緩和時間の両端（τｍｉｎ’およびτｍａｘ’）においてあたかもピークを有するようにし、これに対応する対数緩和時間の最大値および最小値とすればよい。なお、図３では、Ｒkは、対数緩和時間の両端においてピークを有しないが、いずれか一方においてピークを有していれば一方だけを設定し直してもよい。このようにして設定し直した場合、式（Ａ）から得られる自然対数緩和時間分布抵抗Ｒｌは、設定された対数緩和時間の範囲にピークを有することになる。
また、測定したインピーダンススペクトルデータによっては、試料に応じて、いずれか一方にのみピークを有する場合もある。この場合には、そのピークが、測定された対数周波数に相当する対数緩和時間の範囲にあるか否かによって、対数緩和時間の最大値および／または最小値を設定し直すか否かを判定すればよい。
図２および図３を参照して説明したステップＳ１１０の再設定をし直すか否かの判断は、オペレータがディスプレイ等に表示されたデータを見て、目視・手動にて行ってもよいし、コンピュータが判断してもよい。コンピュータが判断し、再設定する場合には、表示されたデータを外挿し、例えば、図２の（Ａ）であれば、インピーダンス実数項のプラトーな対数周波数範囲が、ｆｍａｘ’Ｈｚ以上ｆｍｉｎ’Ｈｚ以下となるようにして、対数緩和時間の最大値および／または最小値を再設定すればよい。
なお、ステップＳ１１０において設定し直すと判定しない場合には、測定された対数周波数に対応する対数緩和時間の最大値および最小値をステップＳ１４０の解析に使用すればよい。
再度図１に戻る。
ステップＳ１４０において、式（Ａ）におけるキャパシタのインピーダンスは、次式（Ｂ）で与えられてもよい。
Ｚ（ｆ）＝Ｔ／（ｊ２πｆ） ・・・（Ｂ）
従来、正則化最小二乗法は、パラメータ値を最小化する項を含むため、キャパシタンス（Ｃ）のインピーダンスを、Ｚ（ｆ）＝１／（Ｃｊ２πｆ）で与えていた。しかしながら、本願発明者は、キャパシタのインピーダンスの式では計算が発散する場合があるが、本願発明者が提案する上記式（Ｂ）を採用することにより、計算の発散を抑制し、正則化非線形最適化を適用できることを見出した。
ステップＳ１４０において、式（Ａ）による非線形最適化は、パラメータ数が多すぎる悪設定問題または不良設定問題になる場合があるため、正則化非線形最適化法を用いて最適化してもよい。このような正則化非線形最適化法は、Ｔｉｋｈｏｎｏｖ正則化法、Ｌａｓｓｏ正則化法およびＥｌａｓｔｉｃ Ｎｅｔ正則化法からなる群から選択される正則化法を採用してもよい。あるいは、正則化非線形最適化法は、Ｌｅｖｅｎｂｅｒｇ－Ｍａｒｑｕａｒｄｔ法、シンプレックス法、および、共役勾配法からなる群から選択される非線形最適化法を採用してもよい。当然ながら、これらの正則化法と非線形最適化法とを組み合わせてもよい。
ステップＳ１４０において、パラメータ初期値は適当な正値を事前に設定するとよい。後述する実施例ではＲ∞＝１、Ｔ＝１０ｆｍａｘ’、Ｌ＝１０ｆｍａｘ’＋３、Ｒｌ＝｜Ｚｉｍａｇ ａｖ｜／Ｍとしている。｜Ｚｉｍａｇ ａｖ｜は測定された全インピーダンス虚数値の平均値の絶対値である。初期値の設定値は０以外の正値かつ正則化非線形最適化の計算によって測定結果を再現しうる値であればよく、上記の与え方に限られるものではない。
ステップＳ１１０～Ｓ１４０によって、パラメータＲ∞、Ｋ、ＬおよびＲｌを高精度に求めることができる。なお、ステップＳ１４０を行った後、再度ステップＳ１２０および／またはステップＳ１３０の設定を見直し、ステップＳ１４０を繰り返し行ってもよい。
本発明によれば、周波数の低下によってインピーダンス絶対値が一定の値に収束しない、または、インピーダンス絶対値が一定の値に収束していないデータに対しても、測定されたインピーダンススペクトルデータを再現する各種パラメータを求めることができるため、測定されたインピーダンススペクトルデータとして、一次電池、二次電池、燃料電池、キャパシタ材料、防蝕処理を施した金属およびセラミックス材料からなる群から選択される試料のインピーダンススペクトルデータを採用できる。
図４は、図１に続く本発明のインピーダンススペクトルデータを解析処理する工程を示すフローチャートである。
ステップＳ１４０：ステップＳ１２０およびステップＳ１３０に続いて、式（Ａ）を満たすようにパラメータを求める。
ステップＳ１５０：ステップＳ１４０で得られたパラメータを式（Ａ）に代入し、理論インピーダンススペクトルデータを求める。
ステップＳ１６０：理論インピーダンススペクトルと測定されたインピーダンススペクトルとを重ねて表示して比較する。
このようにすることにより、ステップＳ１４０で得られたパラメータの妥当性を目視にて瞬時に判断できる。詳細には、ステップＳ１５０で求めた理論インピーダンススペクトルデータと測定されたインピーダンススペクトルデータとが良好に一致すれば、ステップＳ１４０で得られたパラメータが妥当であると判断できる。一方、理論インピーダンススペクトルデータと測定されたインピーダンススペクトルデータとが一致しない場合には、再度、ステップＳ１２０および／またはステップＳ１３０において、対数緩和時間の最大値および／または最小値、ならびに、分割数を設定し直し、ステップＳ１４０、ステップＳ１５０およびステップＳ１６０を繰り返してもよい。
次に、図１に記載の解析処理する方法を実施するインピーダンススペクトルデータ解析システムについて説明する。
図５は、本発明のインピーダンススペクトルデータ解析処理システムを示すブロック図である。
本発明のインピーダンススペクトルデータ解析処理システム（以降では単に解析処理システムと称する）５００は、少なくとも、制御装置５１０と、表示装置５２０とを備える。
制御装置５１０は、数緩和時間判定部５１１において、測定されたインピーダンススペクトルデータに基づいて、測定された対数周波数に対応する対数緩和時間の最大値および／または最小値を設定し直すか否かを判定する。設定し直すと判定した場合は、対数緩和時間設定部５１２において対数緩和時間の最大値および／または最小値を設定し直す。分割数設定部５１３は、設定し直さないと判定した場合には、測定された対数周波数に対応する対数緩和時間の範囲で、または対数緩和時間設定部５１２において設定し直された対数緩和時間の範囲で等間隔に分割する数を設定する（図１のステップＳ１２０およびＳ１３０）。ここで、対数緩和時間の最大値および／または最小値を設定し直すか否かの判定や、設定した対数緩和時間の範囲で等間隔に分割する数の設定の詳細は、図１～図３を参照して説明したとおりであるため省略する。
さらに、制御装置５１０のパラメータ解析部５１４において、正則化最小二乗法を適用し、次式（Ａ）を満たすように、設定した対数緩和時間の範囲および設定した分割数について、パラメータＲ∞、Ｔ、ＬおよびＲｌを解析する（図１のステップＳ１４０）。
Figure 0007026972000010
ここで、ｐは測定されたデータ点の配列番号、ｌは対数緩和時間の配列番号、ｆｐはｐ番目の周波数、ｊは単位虚数、Ｍは対数緩和時間の分割数、Ｔはキャパシタンス（Ｃ）の逆数、Ｌはインダクタンス、Ｒ∞は高周波限抵抗、τｌはｌ番目の緩和時間、Ｒｌはτｌにおける抵抗値、｜Δｌｎτ｜は自然対数緩和時間の間隔の絶対値である。
表示装置５２０は、パラメータ解析部５１４における解析結果を表示する。表示装置５２０は、解析結果を表示できれば任意の表示装置が採用されるが、ディスプレイ装置やプロジェクタ装置である。
制御装置５１０は、入力端末（図示せず）をさらに備えてもよい。これにより、オペレータが、表示装置５２０に表示されたデータを目視にて確認し、測定された対数周波数に対応する対数緩和時間の最大値および／または最小値を設定し直すか否かを判定し、設定した対数緩和時間の範囲で等間隔に分割する数を設定する。
制御装置５１０は、得られたパラメータを、再度式（Ａ）に代入し、理論インピーダンススペクトルデータを求め、理論インピーダンススペクトルデータと、測定されたインピーダンススペクトルデータとを表示装置５２０に表示させ、比較する。
解析処理システム５００は、インピーダンススペクトルデータを測定する電気化学測定装置５３０と接続されていてもよく、測定されたインピーダンススペクトルデータを制御装置５１０が受け取るように構成されてもよい。
また、制御装置５１０は、上述の機能を発揮するものであれば制限はないが、好ましくは、少なくとも、メモリ、中央処理装置（ＣＰＵ）、ハードディスクドライブ（ＨＤＤ）、ドライブ装置、および、通信制御部を備えるコンピュータ装置であり得る。詳細には、オペレーティングシステム（ＯＳ）および図１に示す解析処理を実施するためのアプリケーションプログラム（インピーダンススペクトルデータ解析処理プログラム）がＨＤＤに格納され、ＣＰＵによって実行され、ＨＤＤからメモリに読みだされる。ＣＰＵが、アプリケーションプログラムの内容に応じて、ドライブ装置、通信制御部あるいは表示装置５２０の動作を制御する。
通信制御部が測定されたインピーダンススペクトルデータを受信し、メモリに格納するか、あるいは、電気化学測定装置５３０で測定されたインピーダンススペクトルデータをメモリに格納する。次いで、コンピュータ装置のＣＰＵは、アプリケーションプログラムを実行し、測定されたインピーダンススペクトルデータを読み出し、そのデータに基づいて、測定された対数周波数に対応する対数緩和時間の最大値および／または最小値を設定し直すか否かを判定するステップ（図１のステップＳ１１０）と、判定するステップで設定した対数緩和時間の範囲で等間隔に分割する数を設定するステップ（図１のステップＳ１３０）と、正則化最小二乗法を適用し、次式（Ａ）を満たすように、設定した対数緩和時間の範囲および設定した分割数について、パラメータＲ∞、Ｔ、ＬおよびＲｌを解析するステップ（図１のステップＳ１４０）とを行う。
Figure 0007026972000011
ここで、ｐは測定されたデータ点の配列番号、ｌは対数緩和時間の配列番号、ｆｐはｐ番目の周波数、ｊは単位虚数、Ｍは対数緩和時間の分割数、Ｔはキャパシタンス（Ｃ）の逆数、Ｌはインダクタンス、Ｒ∞は高周波限抵抗、τｌはｌ番目の緩和時間、Ｒｌはτｌにおける抵抗値、｜Δｌｎτ｜は自然対数緩和時間の間隔の絶対値である。
なお、上述のステップＳ１１０～Ｓ１４０は、上述したとおりであるため重複して説明するのを避ける。解析するステップで得られたパラメータは、メモリに格納されるようにしてもよい。
さらに、コンピュータ装置のＣＰＵは、アプリケーションプログラムをさらに実行し、解析するステップで得られたパラメータを式（Ａ）に代入し、理論インピーダンスを求めるステップ（図４のＳ１５０）と、求めた理論インピーダンススペクトルと測定されたインピーダンススペクトルとを表示装置に表示し、比較するステップ（図４のＳ１６０）とを行ってもよい。ここでも、上述のステップＳ１５０～Ｓ１６０は、上述したとおりであるため重複して説明するのを避ける。
上述のアプリケーションプログラム（インピーダンススペクトルデータ解析処理プログラム）は、ＨＤＤに格納されるものとして説明したが、コンピュータ読み取り可能なリムーバルディスク等に格納され頒布される、または、インターネットなどのネットワーク回線や通信制御部を通じて頒布され、コンピュータにインストールされてもよい。
次に具体的な実施例を用いて本発明を詳述するが、本発明がこれら実施例に限定されないことに留意されたい。
［比較例１］
比較例１では、図６に示す電気回路に周波数１ＭＨｚから１Ｈｚ、対数周波数で等間隔に１桁あたり１０点のインピーダンススペクトルデータに基づいて、上述の式（ａ）および（ｂ）を用いて自然対数緩和時間分布抵抗Ｒｋの対数緩和時間依存性（従来法の緩和時間分布解析）を調べた。さらに、得られた緩和時間分布解析から理論インピーダンススペクトルデータを計算し、測定されたインピーダンススペクトルデータと比較した。これらの結果を図７～図９に示す。
図６は、比較例１で用いた電気回路を示す図である。
図７は、図６に示す電気回路を用いて１ＭＨｚから１Ｈｚの周波数で計算したインピーダンススペクトルデータの一覧である。
図７の（Ａ）は、インピーダンス実数項と対数周波数との関係を示し、図７の（Ｂ）は、インピーダンス虚数項と対数周波数との関係を示し、図７の（Ｃ）は、インピーダンスプロットを示す。図７の（Ａ）によれば、インピーダンス実数項は、測定された対数周波数の範囲のうち低周波数側および高周波数側でプラトーでなかった。図７の（Ｂ）によれば、インピーダンス虚数項は、測定された対数周波数の範囲のうち低周波数側および高周波数側でピークを有しなかった。図７の（Ｃ）は、インピーダンスプロットが円弧を有しなかった。
図８は、比較例１による図７のインピーダンススペクトルデータに基づく自然対数緩和時間分布抵抗Ｒｋの対数緩和時間依存性を示す図である。
図８によれば、従来法の緩和時間分布解析では、測定された対数周波数に相当する対数緩和時間の範囲にＲｋのピークを有しないことが分かった。
図９は、比較例１による図８で得られた緩和時間分布解析から計算した理論インピーダンススペクトルデータと、測定されたインピーダンススペクトルデータとの比較を示す図である。
図９の（Ａ）は、インピーダンス実数項と対数周波数との関係を示し、図９の（Ｂ）は、インピーダンス虚数項と対数周波数との関係を示し、図９の（Ｃ）は、インピーダンスプロットを示し、白丸が測定されたインピーダンススペクトルデータ（すなわち図７と同じ）を表し、実線が理論インピーダンススペクトルデータを表す。図９によれば、理論インピーダンススペクトルデータは、とりわけ点線で示す領域において、測定されたインピーダンススペクトルデータを再現できていないことが確認された。
［実施例１］
実施例１では、比較例１と同じ電気回路で測定されたインピーダンススペクトルデータ（図７）に基づいて、図１～図４に記載の解析処理方法を実施した。詳細には、図７の（Ａ）～（Ｃ）、および、図８を参照し、測定された対数周波数に対応する対数緩和時間の最大値および／または最小値を設定し直すか否かを判定した（図１のステップＳ１１０）。判定は目視にて行った。
図７の（Ａ）において、インピーダンス実数項が、測定された対数周波数の範囲（１Ｈｚ～１ＭＨｚ）を超えて低周波数側および高周波数側でプラトーとなるように、図７の（Ｂ）において、インピーダンス虚数項が、測定された対数周波数の範囲（１Ｈｚ～１ＭＨｚ）を超えて低周波数側および高周波数側にピークを有するように、図７の（Ｃ）において、インピーダンスプロットが円弧を有するように、ならびに、図８において、対数緩和時間の範囲を超えてＲｋがあたかもピークを有するように、対数緩和時間の最大値および最小値を、ｌｏｇ（１／（２π×１０－２［ｓ］））およびｌｏｇ（１／（２π×１０９［ｓ］））に手動にて設定し直した。
設定された対数緩和時間の範囲で等間隔に分割する数を対数スケールで緩和時間１桁あたり１０点に手動にて設定した（図１のステップＳ１３０）。次いで、正則化最小二乗法を適用し、式（Ａ）を満たすように、設定した対数緩和時間の範囲および設定した分割数について、パラメータＲ∞、Ｔ、ＬおよびＲｌを解析した（図１のステップＳ１４０）。この解析には、さらに、式（Ｂ）、および、Ｔｉｋｈｏｎｏｖの正則化非線形最適化法を用い、コンピュータ装置が実施した。なお、パラメータの初期値としてＲ∞＝１、Ｔ＝１０ｆｍａｘ’、Ｌ＝１０ｆｍａｘ’＋３、Ｒｋ＝Ｚｉｍａｇ ａｖ／Ｍを用いた。ここで、Ｚｉｍａｇ ａｖは測定された全インピーダンス虚数値の平均値（１．７Ω）であり、ｆｍａｘ’は１０６Ｈｚであり、Ｍは９１であった。結果を図１０に示す。さらに、緩和時間分布解析から理論インピーダンススペクトルデータを計算し、測定されたインピーダンススペクトルデータと比較した。結果を図１１に示す。
図１０は、実施例１による図７のインピーダンススペクトルデータに基づいて、対数緩和時間の範囲を再設定した自然対数緩和時間分布抵抗Ｒｌの対数緩和時間依存性を示す図である。
図１０によれば、対数緩和時間の範囲を再設定し、本発明者が考案した式（Ａ）を用いた緩和時間分布解析では、設定された対数緩和時間の範囲にＲｌのピークを有することが分かった。
図１１は、実施例１による図１０で得られた緩和時間分布解析から計算した理論インピーダンススペクトルデータと、測定されたインピーダンススペクトルデータとの比較を示す図である。
図１１の（Ａ）は、インピーダンス実数項と対数周波数との関係を示し、図１１の（Ｂ）は、インピーダンス虚数項と対数周波数との関係を示し、図１１の（Ｃ）は、インピーダンスプロットを示し、白丸が測定されたインピーダンススペクトルデータ（すなわち図７と同じ）を表し、実線が理論インピーダンススペクトルデータを表す。図１１によれば、理論インピーダンススペクトルデータは、測定されたインピーダンススペクトルデータを完全に再現でき、得られたパラメータが妥当であることが目視にて確認された。
［比較例２］
比較例２では、市販のコイン型リチウム電池（ＣＲ１６２０）の測定されたインピーダンススペクトルデータに基づいて、上述の式（ａ）および（ｂ）を用いて自然対数緩和時間分布抵抗Ｒｋの対数緩和時間依存性（従来法の緩和時間分布解析）を調べた。市販のコイン型リチウム電池（ＣＲ１６２０）のインピーダンススペクトル測定には電気化学測定装置（ソーラートロン製，１２９６型ポテンショ・ガルバノスタットと同社製１２５５型周波数応答特性装置を接続した装置）を用いて測定した。さらに、得られた緩和時間分布解析から理論インピーダンススペクトルデータを計算し、測定されたインピーダンススペクトルデータと比較した。これらの結果を図１２～図１４に示す。
図１２は、比較例２のコイン型リチウム電池のインピーダンススペクトルデータを示す図である。
図１２の（Ａ）は、インピーダンス実数項と対数周波数との関係を示し、図１２の（Ｂ）は、インピーダンス虚数項と対数周波数との関係を示し、図１２の（Ｃ）は、インピーダンスプロットを示すが、周波数が小さくなってもインピーダンス実数項およびインピーダンス虚数項ともに一定値に収束しないスペクトルになった。これは、コイン型リチウム電池において、低周波数域で物質拡散に伴う特性インピーダンスが現れるためである。
図１３は、比較例２による図１２のインピーダンススペクトルデータに基づく自然対数緩和時間分布抵抗Ｒｋの対数緩和時間依存性を示す図である。
図１３によれば、従来法の緩和時間分布解析では、測定された対数周波数に相当する対数緩和時間の範囲にＲｋのピークを有しないことが分かった。
図１４は、比較例２による図１３で得られた緩和時間分布解析から計算した理論インピーダンススペクトルデータと、測定されたインピーダンススペクトルデータとの比較を示す図である。
図１４の（Ａ）は、インピーダンス実数項と対数周波数との関係を示し、図１４の（Ｂ）は、インピーダンス虚数項と対数周波数との関係を示し、図１４の（Ｃ）は、インピーダンスプロットを示し、白丸が測定されたインピーダンススペクトルデータ（すなわち図１２と同じ）を表し、実線が理論インピーダンススペクトルデータを表す。図１４によれば、理論インピーダンススペクトルデータは、とりわけ点線で示す領域（対数周波数ｌｏｇｆが約－１）において、測定されたインピーダンススペクトルデータを再現できていないことが確認された。
［実施例２］
実施例２では、比較例２と同じコイン型リチウム電池で測定されたインピーダンススペクトルデータ（図１２）に基づいて、図１～図４に記載の解析処理方法を実施した。詳細には、図１３を参照し、測定された対数周波数に対応する対数緩和時間の最大値および／または最小値を設定し直すか否かを判定した（図１のステップＳ１１０）。
図１３において、対数緩和時間の範囲を超えてＲｋがあたかもピークを有するように、対数緩和時間の最大値および最小値を、ｌｏｇ（１／（２π×１０－３［ｓ］））およびｌｏｇ（１／（２π×１０６［ｓ］））に手動にて設定し直した。
設定された対数緩和時間の範囲で等間隔に分割する数を対数スケールで緩和時間１桁あたり１０点に手動にて設定した（図１のステップＳ１３０）。次いで、正則化最小二乗法を適用し、式（Ａ）を満たすように、設定した対数緩和時間の範囲および設定した分割数について、パラメータＲ∞、Ｔ、ＬおよびＲｌを解析した（図１のステップＳ１４０）。なお、パラメータの初期値としてＲ∞＝１、Ｔ＝１０ｆｍａｘ’、Ｌ＝１０ｆｍａｘ’＋３、Ｒｋ＝｜Ｚｉｍａｇ ａｖ｜／Ｍを用いた。ここで、｜Ｚｉｍａｇ ａｖ｜は測定された全インピーダンス虚数値の平均値の絶対値（１．７Ω）であり、ｆｍａｘ’は１０６Ｈｚであり、Ｍは９１であった。この解析には、式（Ｂ）、および、Ｔｉｋｈｏｎｏｖの正則化非線形最適化法を用い、コンピュータ装置が実施した。結果を図１５に示す。さらに、緩和時間分布解析から理論インピーダンススペクトルデータを計算し、測定されたインピーダンススペクトルデータと比較した。結果を図１６に示す。
図１５は、実施例２による図１２のインピーダンススペクトルデータに基づいて、対数緩和時間の範囲を再設定した自然対数緩和時間分布抵抗Ｒｌの対数緩和時間依存性を示す図である。
図１５によれば、対数緩和時間の範囲を再設定し、本発明者が考案した式（Ａ）を用いた緩和時間分布解析では、設定された対数緩和時間の範囲にＲｌのピークを有することが分かった。
図１６は、実施例２による図１５で得られた緩和時間分布解析から計算した理論インピーダンススペクトルデータと、測定されたインピーダンススペクトルデータとの比較を示す図である。
図１６の（Ａ）は、インピーダンス実数項と対数周波数との関係を示し、図１６の（Ｂ）は、インピーダンス虚数項と対数周波数との関係を示し、図１６の（Ｃ）は、インピーダンスプロットを示し、白丸が測定されたインピーダンススペクトルデータ（すなわち図１２と同じ）を表し、実線が理論インピーダンススペクトルデータを表す。図１６によれば、理論インピーダンススペクトルデータは、測定されたインピーダンススペクトルデータを完全に再現でき、得られたパラメータが妥当であることが目視にて確認された。
［比較例３］
比較例３では、Ｇｏｏｇｌｅ Ｉｍｐｅｄａｎｃｅ Ｓｐｅｃｔｒｏｓｃｏｐｙフォーラム（ｈｔｔｐｓ：／／ｇｒｏｕｐｓ．ｇｏｏｇｌｅ．ｃｏｍ／ｆｏｒｕｍ／＃！ｆｏｒｕｍ／ｉｍｐｅｄａｎｃｅ－ｓｐｅｃｔｒｏｓｃｏｐｙ）にアップロードされている絶縁性酸化物のインピーダンススペクトル（図１７）に基づいて、上述の式（ａ）および（ｂ）を用いて自然対数緩和時間分布抵抗Ｒｋの対数緩和時間依存性（従来法の緩和時間分布解析）を調べた。さらに、得られた緩和時間分布解析から理論インピーダンススペクトルデータを計算し、測定されたインピーダンススペクトルデータと比較した。これらの結果を図１８および図１９に示す。
図１７は、比較例３のＧｏｏｇｌｅ Ｉｍｐｅｄａｎｃｅ Ｓｐｅｃｔｒｏｓｃｏｐｙフォーラムにアップロードされた絶縁性酸化物のインピーダンススペクトルデータを示す図である。
図１７の（Ａ）は、インピーダンス実数項と対数周波数との関係を示し、図１７の（Ｂ）は、インピーダンス虚数項と対数周波数との関係を示し、図１７の（Ｃ）は、インピーダンスプロットを示すが、周波数が小さくなってもインピーダンス実数項およびインピーダンス虚数項ともに一定値に収束しないスペクトルになった。
図１８は、比較例３による図１７のインピーダンススペクトルデータに基づく自然対数緩和時間分布抵抗Ｒｋの対数緩和時間依存性を示す図である。
図１８によれば、従来法の緩和時間分布解析では、測定された対数周波数に相当する対数緩和時間の範囲にＲｋのピークを有しないことが分かった。
図１９は、比較例３による図１８で得られた緩和時間分布解析から計算した理論インピーダンススペクトルデータと、アップロードされたインピーダンススペクトルデータとの比較を示す図である。
図１９の（Ａ）は、インピーダンス実数項と対数周波数との関係を示し、図１９の（Ｂ）は、インピーダンス虚数項と対数周波数との関係を示し、図１９の（Ｃ）は、インピーダンスプロットを示し、白丸がアップロードされたインピーダンススペクトルデータ（すなわち図１７と同じ）を表し、実線が理論インピーダンススペクトルデータを表す。図１９によれば、理論インピーダンススペクトルデータは、とりわけ点線で示す領域において、測定されたインピーダンススペクトルデータを再現できていないことが確認された。
［実施例３］
実施例３では、比較例３と用いたインピーダンススペクトルデータ（図１７）に基づいて、図１～図４に記載の解析処理方法を実施した。詳細には、図１８を参照し、測定された対数周波数に対応する対数緩和時間の最大値および／または最小値を設定し直すか否かを判定した（図１のステップＳ１１０）。
図１８において、対数緩和時間の範囲を超えてＲｋがあたかもピークを有するように、対数緩和時間の最大値および最小値を、ｌｏｇ（１／（２π×１０－２［ｓ］））およびｌｏｇ（１／（２π×１０６［ｓ］））に手動にて設定し直した。
設定された対数緩和時間の範囲で等間隔に分割する数を対数スケールで緩和時間１桁あたり１０点に手動にて設定した（図１のステップＳ１３０）。次いで、正則化最小二乗法を適用し、式（Ａ）を満たすように、設定した対数緩和時間の範囲および設定した分割数について、パラメータＲ∞、Ｔ、ＬおよびＲｌを解析した（図１のステップＳ１４０）。なお、パラメータの初期値としてＲ∞＝１、Ｔ＝１０ｆｍａｘ’、Ｌ＝１０ｆｍａｘ’＋３、Ｒｋ＝｜Ｚｉｍａｇ ａｖ｜／Ｍを用いた。ここで、｜Ｚｉｍａｇ ａｖ｜は測定された全インピーダンス虚数値の平均値の絶対値（１．２×１０９Ω）であり、ｆｍａｘ’は１０６Ｈｚであり、Ｍは７１であった。この解析には、式（Ｂ）、および、Ｔｉｋｈｏｎｏｖの正則化非線形最適化法を用い、コンピュータ装置が実施した。結果を図２０に示す。さらに、緩和時間分布解析から理論インピーダンススペクトルデータを計算し、アップロードされたインピーダンススペクトルデータと比較した。結果を図２１に示す。
図２０は、実施例３による図１７のインピーダンススペクトルデータに基づいて、対数緩和時間の範囲を再設定した自然対数緩和時間分布抵抗Ｒｌの対数緩和時間依存性を示す図である。
図２０によれば、対数緩和時間の範囲を再設定し、本発明者が考案した式（Ａ）を用いた緩和時間分布解析では、設定された対数緩和時間の範囲にＲｌのピークを有することが分かった。
図２１は、実施例３による図２０で得られた緩和時間分布解析から計算した理論インピーダンススペクトルデータと、アップロードされたインピーダンススペクトルデータとの比較を示す図である。
図２１の（Ａ）は、インピーダンス実数項と対数周波数との関係を示し、図２１の（Ｂ）は、インピーダンス虚数項と対数周波数との関係を示し、図２１の（Ｃ）は、インピーダンスプロットを示し、白丸がアップロードされたインピーダンススペクトルデータ（すなわち図１７同じ）を表し、実線が理論インピーダンススペクトルデータを表す。図２１によれば、理論インピーダンススペクトルデータは、アップロードされたインピーダンススペクトルデータを完全に再現でき、得られたパラメータが妥当であることが目視にて確認された。
本発明のインピーダンススペクトルデータ解析処理方法を採用すれば、インピーダンス絶対値が一定の値に収束しないようなデータに対してもモデルによる再現を可能とし、緩和時間分布解析も可能となる。このような解析処理方法を採用すれば、インピーダンススペクトルデータ解析処理システムおよび解析処理プログラムを提供でき、一次電池、二次電池、キャパシタ、腐食系など種々の試料に対してインピーダンススペクトルデータの解析が可能となる。
５００ インピーダンススペクトルデータ解析処理システム
５１０ 制御装置
５２０ 表示装置
５３０ 電気化学測定装置

Claims (21)

  1. 測定されたインピーダンススペクトルデータを用いて解析処理する方法であって、
    前記測定されたインピーダンススペクトルデータに基づいて、測定された対数周波数に対応する対数緩和時間の最大値および／または最小値を設定し直すか否かを判定するステップと、
    前記判定するステップで設定し直すと判定した場合、前記対数緩和時間の最大値および／または最小値を設定し直すステップと、
    前記判定するステップで設定し直さないと判定した場合には前記測定された対数周波数に対応する対数緩和時間の範囲、または、前記判定するステップで設定し直すと判定した場合には前記設定し直された対数緩和時間の範囲で等間隔に分割する分割数を設定するステップと、
    正則化最小二乗法を適用し、次式（Ａ）
    Figure 0007026972000012

    （ｐは測定されたデータ点の配列番号、ｌは対数緩和時間の配列番号、ｆｐはｐ番目の周波数、ｊは単位虚数、Ｍは対数緩和時間の分割数、Ｔはキャパシタンス（Ｃ）の逆数、Ｌはインダクタンス、Ｒ∞は高周波限抵抗、τｌはｌ番目の緩和時間、Ｒｌはτｌにおける抵抗値、｜Δｌｎτ｜は自然対数緩和時間の間隔の絶対値）
    を満たすように、前記対数緩和時間の範囲および設定した前記分割数について、パラメータＲ∞、Ｔ、ＬおよびＲｌを解析するステップと
    を包含する方法。
  2. 前記判定するステップは、前記測定されたインピーダンススペクトルデータが、
    （１）インピーダンス実数項と対数周波数との関係において、インピーダンス実数項が、測定された対数周波数の範囲のうち低周波数側および高周波数側でプラトーを有さないか、
    （２）インピーダンス虚数項と対数周波数との関係において、インピーダンス虚数項が、測定された対数周波数の範囲のうち低周波数側および高周波数側でピークを有さないか、または、
    （３）インピーダンスプロットにおいて円弧を有さない場合に
    前記解析するステップに用いる対数緩和時間の最大値および／または最小値を設定し直すと判定し、
    前記設定し直すステップは、
    （１’）前記インピーダンス実数項と対数周波数との関係において、前記インピーダンス実数項が、前記測定された対数周波数の範囲を超えて低周波数側および高周波数側でプラトーとなるように前記解析するステップに用いる対数緩和時間の最大値および／または最小値を設定し直すか、
    （２’）前記インピーダンス虚数項と対数周波数との関係において、前記インピーダンス虚数項が、前記測定された対数周波数の範囲を超えて低周波数側および高周波数側でピークを有するように前記解析するステップに用いる対数緩和時間の最大値および／または最小値を設定し直すか、
    （３’）前記インピーダンスプロットにおいて円弧を有するように前記解析するステップに用いる対数緩和時間の最大値および／または最小値を設定し直す
    請求項１に記載の方法。
  3. 前記判定するステップは、前記測定された対数周波数の範囲について、次式（ａ）および次式（ｂ）
    Figure 0007026972000013

    Figure 0007026972000014

    （ｋは対数緩和時間の配列番号、ｐは測定されたデータ点の配列番号、ｆｋおよびｆｐはｋ番目およびｐ番目の周波数、ｊは単位虚数、Ｎは測定データ点数、Ｒ∞は高周波限抵抗、τｋはｋ番目の緩和時間、Ｒｋはτｋにおける抵抗値、｜Δｌｎτ｜は自然対数緩和時間の間隔の絶対値）
    を用いて、自然対数緩和時間分布抵抗Ｒｋの対数緩和時間依存性を求め、前記測定された対数周波数に相当する対数緩和時間の範囲に前記Ｒｋがピークを有さない場合に、前記解析するステップに用いる対数緩和時間の最大値および／または最小値を設定し直すと判定し、
    前記設定し直すステップは、前記測定された対数周波数に相当する対数緩和時間の範囲を超えて前記Ｒｋがピークを有するように、前記解析するステップに用いる対数緩和時間の最大値および／または最小値を設定し直す
    請求項１に記載の方法。
  4. 前記解析するステップにおいて、前記式（Ａ）における前記キャパシタンスのインピーダンスは、次式（Ｂ）
    Ｚ（ｆ）＝Ｔ／（ｊ２πｆ） ・・・（Ｂ）
    で与えられる、請求項１～３のいずれかに記載の方法。
  5. 前記解析するステップは、正則化非線形最適化法を用いて算出する、請求項１～４のいずれかに記載の方法。
  6. 前記正則化非線形最適化法は、Ｔｉｋｈｏｎｏｖ正則化法、Ｌａｓｓｏ正則化法およびＥｌａｓｔｉｃ Ｎｅｔ正則化法からなる群から選択される正則化法を採用する、請求項５に記載の方法。
  7. 前記正則化非線形最適化法は、Ｌｅｖｅｎｂｅｒｇ－Ｍａｒｑｕａｒｄｔ法、シンプレックス法、および、共役勾配法からなる群から選択される非線形最適化法を採用する、請求項５または６に記載の方法。
  8. 前記解析するステップで得られた前記パラメータを前記式（Ａ）に代入し、理論インピーダンススペクトルを求めるステップと、
    前記求めた理論インピーダンススペクトルと前記測定されたインピーダンススペクトルとを表示し、比較し、一致するか否かを判定するステップと
    をさらに包含する、請求項１～７のいずれかに記載の方法。
  9. 前記表示し、比較し、一致するか否かを判定するステップが、一致しないと判定した場合、前記設定し直すステップ以降を繰り返す、請求項８に記載の方法。
  10. 前記測定されたインピーダンススペクトルデータは、一次電池、二次電池、燃料電池、キャパシタ材料、防蝕処理を施した金属およびセラミックス材料からなる群から選択される試料のインピーダンススペクトルデータである、請求項１～９のいずれかに記載の方法。
  11. 測定されたインピーダンススペクトルデータを用いたインピーダンススペクトルデータ解析処理システムであって、
    制御装置と、
    表示装置と
    を備え、
    前記制御装置は、
    前記測定されたインピーダンススペクトルデータに基づいて、測定された対数周波数に対応する対数緩和時間の最大値および／または最小値を設定し直すか否かを判定する対数緩和時間判定部と、
    前記対数緩和時間判定部が設定し直すと判定した場合、前記測定された対数周波数に対応する対数緩和時間の最大値および／または最小値を設定し直す対数緩和時間設定部と、
    前記対数緩和時間判定部が設定し直さないと判定した場合には前記測定された対数周波数に対応する対数緩和時間の範囲、または、前記対数緩和時間判定部が設定し直すと判定し場合には前記対数緩和時間設定部が設定し直した対数緩和時間の範囲で等間隔に分割する分割数を設定する分割数設定部と、
    正則化最小二乗法を適用し、次式（Ａ）
    Figure 0007026972000015

    （ｐは測定されたデータ点の配列番号、ｌは対数緩和時間の配列番号、ｆｐはｐ番目の周波数、ｊは単位虚数、Ｍは対数緩和時間の分割数、Ｔはキャパシタンス（Ｃ）の逆数、Ｌはインダクタンス、Ｒ∞は高周波限抵抗、τｌはｌ番目の緩和時間、Ｒｌはτｌにおける抵抗値、｜Δｌｎτ｜は自然対数緩和時間の間隔の絶対値）
    を満たすように、前記対数緩和時間の範囲および設定した前記分割数について、パラメータＲ∞、Ｔ、ＬおよびＲｌを解析するパラメータ解析部と
    をさらに備え、
    前記表示装置は、前記制御装置における解析結果を表示する
    解析処理システム。
  12. 前記制御装置の前記対数緩和時間判定部は、前記測定されたインピーダンススペクトルデータが、
    （１）インピーダンス実数項と対数周波数との関係において、インピーダンス実数項が、測定された対数周波数の範囲のうち低周波数側および高周波数側でプラトーを有さないか、
    （２）インピーダンス虚数項と対数周波数との関係において、インピーダンス虚数項が、測定された対数周波数の範囲のうち低周波数側および高周波数側でピークを有さないか、または、
    （３）インピーダンスプロットにおいて円弧を有さない場合に
    前記パラメータ解析部で用いる対数緩和時間の最大値および／または最小値を設定し直すと判定し、
    前記制御装置の前記対数緩和時間設定部は、
    （１’）前記インピーダンス実数項と対数周波数との関係において、前記インピーダンス実数項が、前記測定された対数周波数の範囲を超えて低周波数側および高周波数側でプラトーとなるように、前記パラメータ解析部で用いる対数緩和時間の最大値および／または最小値を設定し直すか、
    （２’）前記インピーダンス虚数項と対数周波数との関係において、前記インピーダンス虚数項が、前記測定された対数周波数の範囲を超えて低周波数側および高周波数側でピークを有するように、前記パラメータ解析部で用いる対数緩和時間の最大値および／または最小値を設定し直すか、
    （３’）前記インピーダンスプロットにおいて円弧を有するように前記パラメータ解析部で用いる対数緩和時間の最大値および／または最小値を設定し直す
    請求項１１に記載の解析処理システム。
  13. 前記制御装置の前記対数緩和時間判定部は、前記測定された対数周波数の範囲について、次式（ａ）および次式（ｂ）
    Figure 0007026972000016

    Figure 0007026972000017

    （ｋは対数緩和時間の配列番号、ｐは測定されたデータ点の配列番号、ｆｋおよびｆｐはｋ番目およびｐ番目の周波数、ｊは単位虚数、Ｎは測定データ点数、Ｒ∞は高周波限抵抗、τｋはｋ番目の緩和時間、Ｒｋはτｋにおける抵抗値、｜Δｌｎτ｜は自然対数緩和時間の間隔の絶対値）
    を用いて、自然対数緩和時間分布抵抗Ｒｋの対数緩和時間依存性を求め、前記測定された対数周波数の範囲に前記Ｒｋがピークを有さない場合に、前記パラメータ解析部で用いる対数緩和時間の最大値および／または最小値を設定し直すと判定し、
    前記制御装置の前記対数緩和時間設定部は、前記測定された対数周波数に相当する対数緩和時間の範囲を超えて前記Ｒｋがピークを有するように前記パラメータ解析部で用いる対数緩和時間の最大値および／または最小値を設定し直す
    請求項１１に記載の解析処理システム。
  14. 前記制御装置は、オペレータが、前記表示装置に表示された前記測定されたインピーダンススペクトルデータから目視にて前記（１）～（３）を判定し、前記（１’）～（３’）を手動にて設定し直す入力端末を備える、請求項１２に記載の解析処理システム。
  15. 前記制御装置は、オペレータが、前記表示装置に表示された前記自然対数緩和時間分布抵抗Ｒｋの対数緩和時間依存性から目視にて前記測定された対数周波数の範囲に前記Ｒｋがピークを有するか否かを判定し、前記測定された対数周波数の範囲を超えて前記Ｒｋがピークを有するように、前記パラメータ解析部で用いる対数緩和時間の最大値および／または最小値を手動にて設定し直す入力端末を備える、請求項１３に記載の解析処理システム。
  16. 前記式（Ａ）において、前記キャパシタンスのインピーダンスは、次式（Ｂ）
    Ｚ（ｆ）＝Ｔ／（ｊ２πｆ） ・・・（Ｂ）
    で与えられる、請求項１１～１５のいずれかに記載の解析処理システム。
  17. 前記制御装置は、得られた前記パラメータを前記式（Ａ）に代入し、理論インピーダンススペクトルを求め、前記求めた理論インピーダンススペクトルと前記測定されたインピーダンススペクトルとを前記表示装置に表示し、比較し、前記理論インピーダンススペクトルと前記測定されたインピーダンススペクトルとが一致するか否かを判定する判定部をさらに備える、請求項１１～１６のいずれかに記載の解析処理システム。
  18. コンピュータに、
    測定されたインピーダンススペクトルデータに基づいて、測定された対数周波数に対応する対数緩和時間の最大値および／または最小値を設定し直すか否かを判定するステップと、
    前記判定するステップで設定し直すと判定した場合、前記対数緩和時間の最大値および／または最小値を設定し直すステップと、
    前記判定するステップで設定し直さないと判定した場合には前記測定された対数周波数に対応する対数緩和時間の範囲、または、前記判定するステップで設定し直すと判定した場合には前記設定し直された対数緩和時間の範囲で等間隔に分割する数を設定するステップと、
    正則化最小二乗法を適用し、次式（Ａ）
    Figure 0007026972000018

    （ｐは測定されたデータ点の配列番号、ｌは対数緩和時間の配列番号、ｆｐはｐ番目の周波数、ｊは単位虚数、Ｍは対数緩和時間の分割数、Ｔはキャパシタンス（Ｃ）の逆数、Ｌはインダクタンス、Ｒ∞は高周波限抵抗、τｌはｌ番目の緩和時間、Ｒｌはτｌにおける抵抗値、｜Δｌｎτ｜は自然対数緩和時間の間隔の絶対値）
    を満たすように、前記対数緩和時間の範囲および設定した分割数について、パラメータＲ∞、Ｔ、ＬおよびＲｌを解析するステップと
    を実行させるための、インピーダンススペクトルデータ解析処理プログラム。
  19. 前記判定するステップは、前記測定されたインピーダンススペクトルデータが、
    （１）インピーダンス実数項と対数周波数との関係において、インピーダンス実数項が、測定された対数周波数の範囲のうち低周波数側および高周波数側でプラトーを有さないか、
    （２）インピーダンス虚数項と対数周波数との関係において、インピーダンス虚数項が、測定された対数周波数の範囲のうち低周波数側および高周波数側でピークを有さないか、または、
    （３）インピーダンスプロットにおいて円弧を有さない場合に、
    前記解析するステップに用いる対数緩和時間の最大値および／または最小値を設定し直すと判定し、
    前記設定し直すステップは、
    （１’）前記インピーダンス実数項と対数周波数との関係において、前記インピーダンス実数項が、前記測定された対数周波数の範囲を超えて低周波数側および高周波数側でプラトーとなるように前記解析するステップに用いる対数緩和時間の最大値および／または最小値を設定し直すか、
    （２’）前記インピーダンス虚数項と対数周波数との関係において、前記インピーダンス虚数項が、前記測定された対数周波数の範囲を超えて低周波数側および高周波数側でピークを有するように前記解析するステップに用いる対数緩和時間の最大値および／または最小値を設定し直すか、
    （３’）前記インピーダンスプロットにおいて円弧を有するように前記解析するステップに用いる対数緩和時間の最大値および／または最小値を設定し直す
    請求項１８に記載のインピーダンススペクトルデータ解析処理プログラム。
  20. 前記判定するステップは、測定された対数周波数の範囲について、次式（ａ）および次式（ｂ）
    Figure 0007026972000019

    Figure 0007026972000020

    （ｋは対数緩和時間の配列番号、ｐは測定されたデータ点の配列番号、ｆｋおよびｆｐはｋ番目およびｐ番目の周波数、ｊは単位虚数、Ｎは測定データ点数、Ｒ∞は高周波限抵抗、τｋはｋ番目の緩和時間、Ｒｋはτｋにおける抵抗値、｜Δｌｎτ｜は自然対数緩和時間の間隔の絶対値）
    を用いて、自然対数緩和時間分布抵抗Ｒｋの対数緩和時間依存性を求め、前記測定された対数周波数の範囲に前記Ｒｋがピークを有さない場合に、前記解析するステップに用いる対数緩和時間の最大値および／または最小値を設定し直すと判定し、
    前記設定し直すステップは、前記測定された対数周波数に相当する対数緩和時間の範囲を超えて前記Ｒｋがピークを有するように前記解析するステップに用いる対数緩和時間の最大値および／または最小値を設定し直す
    請求項１８に記載のインピーダンススペクトルデータ解析処理プログラム。
  21. 前記解析するステップで得られた前記パラメータを前記式（Ａ）に代入し、理論インピーダンススペクトルを求めるステップと、
    前記求めた理論インピーダンススペクトルと前記測定されたインピーダンススペクトルとを表示装置に表示し、比較し、一致するか否かを判定するステップと
    をさらに実行させる、請求項１８～２０のいずれかに記載のインピーダンススペクトルデータ解析処理プログラム。
JP2020532307A 2018-07-23 2019-07-16 インピーダンススペクトルデータを用いた解析処理方法、インピーダンススペクトルデータ解析処理システム、および、インピーダンススペクトル解析処理プログラム Active JP7026972B2 (ja)

Applications Claiming Priority (3)

Application Number Priority Date Filing Date Title
JP2018137404 2018-07-23
JP2018137404 2018-07-23
PCT/JP2019/027936 WO2020022124A1 (ja) 2018-07-23 2019-07-16 インピーダンススペクトルデータを用いた解析処理方法、インピーダンススペクトルデータ解析処理システム、および、インピーダンススペクトル解析処理プログラム

Publications (2)

Publication Number Publication Date
JPWO2020022124A1 JPWO2020022124A1 (ja) 2021-08-19
JP7026972B2 true JP7026972B2 (ja) 2022-03-01

Family

ID=69181418

Family Applications (1)

Application Number Title Priority Date Filing Date
JP2020532307A Active JP7026972B2 (ja) 2018-07-23 2019-07-16 インピーダンススペクトルデータを用いた解析処理方法、インピーダンススペクトルデータ解析処理システム、および、インピーダンススペクトル解析処理プログラム

Country Status (4)

Country Link
US (1) US20210293867A1 (ja)
EP (1) EP3848713B1 (ja)
JP (1) JP7026972B2 (ja)
WO (1) WO2020022124A1 (ja)

Families Citing this family (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN112540316B (zh) * 2020-11-03 2022-07-05 华中科技大学 一种复杂电池阻抗谱解析方法
CN112816895B (zh) * 2020-12-31 2024-04-12 上海派能能源科技股份有限公司 电化学阻抗谱的分析方法、系统、设备及计算机存储介质
CN113777410A (zh) * 2021-09-02 2021-12-10 吉林大学 一种测量BaTiO3介电常数的方法
CN114217136B (zh) * 2022-02-22 2022-05-06 山东卓朗检测股份有限公司 基于大数据的防雷接地电阻检测统计方法
CN115982914B (zh) * 2023-02-16 2023-06-02 华东交通大学 基于弹性网理论的接地网腐蚀评估方法及系统
CN116819377B (zh) * 2023-06-28 2024-05-07 上海炙云新能源科技有限公司 电池组异常识别诊断方法和系统

Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN104914312A (zh) 2015-06-18 2015-09-16 哈尔滨工业大学 一种计算交流阻抗谱弛豫时间分布的方法
WO2016071801A1 (en) 2014-11-04 2016-05-12 Universita' Degli Studi Di Salerno Method and apparatus for monitoring and diagnosing electrochemical devices based on automatic electrochemical impedance identification
WO2018131071A1 (ja) 2017-01-10 2018-07-19 日産自動車株式会社 燃料電池システムの制御方法及び燃料電池システム

Family Cites Families (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US5486762A (en) * 1992-11-02 1996-01-23 Schlumberger Technology Corp. Apparatus including multi-wait time pulsed NMR logging method for determining accurate T2-distributions and accurate T1/T2 ratios and generating a more accurate output record using the updated T2-distributions and T1/T2 ratios
US6321119B1 (en) * 1997-09-24 2001-11-20 Healthonics, Inc. Pulsed signal generator for bioelectric stimulation and healing acceleration
JP3885356B2 (ja) * 1998-04-20 2007-02-21 大日本インキ化学工業株式会社 制振材料の試験方法
CN102544496A (zh) * 2005-09-09 2012-07-04 A123系统公司 具有高充电和放电倍率能力和低阻抗增长的锂二次电池

Patent Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
WO2016071801A1 (en) 2014-11-04 2016-05-12 Universita' Degli Studi Di Salerno Method and apparatus for monitoring and diagnosing electrochemical devices based on automatic electrochemical impedance identification
CN104914312A (zh) 2015-06-18 2015-09-16 哈尔滨工业大学 一种计算交流阻抗谱弛豫时间分布的方法
WO2018131071A1 (ja) 2017-01-10 2018-07-19 日産自動車株式会社 燃料電池システムの制御方法及び燃料電池システム

Non-Patent Citations (2)

* Cited by examiner, † Cited by third party
Title
小林 清（KOBAYASHI,Kiyoshi），緩和時間分布（ＤＲＴ）を用いたＳＯＦＣ電極反応解析の不確定性 Uncertainty of SOFC electrode rcaction analysis by the DRT technique，電気化学会講演要旨集 第８１回大会 ，日本，公益社団法人電気化学会，2014年05月12日，Ｐ。１０６
鷲見 裕史，平野 隆行，緩和時間分布法による等価回路作成支援ソフトウェアの開発，電池討論会講演要旨集 THE 57TH BATTERY SYMPOSIUM IN JAPAN ，日本，（公社）電気化学会 電池技術委員会 第５７回 電池討論会実行委員長 金村 聖志，2016年11月28日，p.195

Also Published As

Publication number Publication date
JPWO2020022124A1 (ja) 2021-08-19
EP3848713A1 (en) 2021-07-14
WO2020022124A1 (ja) 2020-01-30
EP3848713A4 (en) 2022-04-13
US20210293867A1 (en) 2021-09-23
EP3848713B1 (en) 2022-08-24

Similar Documents

Publication Publication Date Title
JP7026972B2 (ja) インピーダンススペクトルデータを用いた解析処理方法、インピーダンススペクトルデータ解析処理システム、および、インピーダンススペクトル解析処理プログラム
Alavi et al. Time-domain fitting of battery electrochemical impedance models
Sun et al. Application of the positive reality principle to metal electrode linear polarization phenomena
JP2000009817A (ja) インピダンススペクトル分析による電池容量測定方法及びその装置
US9194900B2 (en) Electrode evaluation apparatus and electrode evaluation method
Astafev et al. Measurement of electrochemical noise of a Li/MnO 2 primary lithium battery
JP5068383B1 (ja) 予測装置、予測方法及び予測プログラム
DE112014004166T5 (de) Gerät und Verfahren für die Abschätzung der Lebensdauer einer Energiespeichervorrichtung und eines Energiespeichersystems
US11287481B2 (en) Methods and algorithms of cyclic coulometry
TWI759027B (zh) 電池性能評鑑方法及電池性能評鑑裝置
CN112098756A (zh) 电磁兼容问题定位方法、装置、设备及存储介质
CN109116108B (zh) 一种频谱密度图的显示装置和方法
CN111044905A (zh) 数据确定方法及装置
US7577536B1 (en) Determination of ohmic losses in electrical devices
KR20210149795A (ko) 이차전지 검사 방법 및 이차전지 검사 장치
JP6770739B2 (ja) 等価回路推定方法
CN110263417A (zh) 一种时序特性获取方法、装置及电子设备
US20210181290A1 (en) Method of Calibrating Impedance Measurements of A Battery
US10436873B1 (en) Method of calibrating impedance measurements of a battery
EP3980796A1 (en) Systems apparatus and methods of cyclic coulometry
US20130108163A1 (en) Image evaluation apparatus, image evaluation method, and program
CN114711638A (zh) 烘焙设备的控制方法、装置、烘焙设备及存储介质
CN106169944B (zh) 参数自动调整系统及其参数自动调整的方法
CN114062441A (zh) 状态测定装置、状态测定方法及记录介质
WO2019113750A1 (zh) 一种配置用户设备测量参数的方法、用户设备及网络设备

Legal Events

Date Code Title Description
A621 Written request for application examination

Free format text: JAPANESE INTERMEDIATE CODE: A621

Effective date: 20210107

A131 Notification of reasons for refusal

Free format text: JAPANESE INTERMEDIATE CODE: A131

Effective date: 20220107

A521 Request for written amendment filed

Free format text: JAPANESE INTERMEDIATE CODE: A523

Effective date: 20220120

TRDD Decision of grant or rejection written
A01 Written decision to grant a patent or to grant a registration (utility model)

Free format text: JAPANESE INTERMEDIATE CODE: A01

Effective date: 20220202

A61 First payment of annual fees (during grant procedure)

Free format text: JAPANESE INTERMEDIATE CODE: A61

Effective date: 20220208

R150 Certificate of patent or registration of utility model

Ref document number: 7026972

Country of ref document: JP

Free format text: JAPANESE INTERMEDIATE CODE: R150