2023-90638 モータの振動現象の物理量の算出方法、モータ制御方法、算出装置、およびプログラム

(19)【発行国】日本国特許庁(JP)

(12)【公報種別】公開特許公報(A)

(11)【公開番号】P2023090638

(43)【公開日】2023-06-29

(54)【発明の名称】モータの振動現象の物理量の算出方法、モータ制御方法、算出装置、およびプログラム

(51)【国際特許分類】

   G01R 31/34 20200101AFI20230622BHJP

   H02P 23/14 20060101ALI20230622BHJP

【ＦＩ】

G01R31/34 G

H02P23/14

【審査請求】未請求

【請求項の数】10

【出願形態】ＯＬ

(21)【出願番号】P 2022174711

(22)【出願日】2022-10-31

(31)【優先権主張番号】P 2021205483

(32)【優先日】2021-12-17

(33)【優先権主張国・地域又は機関】JP

【新規性喪失の例外の表示】特許法第３０条第２項適用申請有り ２０２２年８月３日に日本能率協会主催テクノフロンティア 第４２回 モータ技術シンポジウムの発表用資料 〔刊行物等〕 ２０２２年８月２３日に一般社団法人電気学会が発行した論文集「２０２２年電気学会産業応用部門大会予稿集」にて掲載

(71)【出願人】

【識別番号】000004204

【氏名又は名称】日本精工株式会社

(74)【代理人】

【識別番号】110002000

【氏名又は名称】弁理士法人栄光事務所

(72)【発明者】

【氏名】新田 勇

【テーマコード（参考）】

2G116

5H505

【Ｆターム（参考）】

2G116BA03

2G116BB06

2G116BB08

5H505BB04

5H505DD08

5H505HB01

5H505JJ03

5H505JJ17

5H505LL01

5H505LL22

5H505LL40

5H505LL41

5H505LL50

(57)【要約】

【課題】従来よりもさらに計算負荷を抑制し、計算の精度を向上させた、モータにおける振動現象の物理量の算出方法を提供する。
【解決手段】モータの振動現象の物理量を算出する算出方法であって、前記モータが備えるステータコアの突極に作用する時間依存の加振力を、空間および時間の周波数ごとに変換して前記モータの振動現象の物理量を算出する算出工程を有し、前記算出工程において、前記ステータコアの突極数に分割した領域のそれぞれに作用する値を前記加振力として用いる。
【選択図】図９

【特許請求の範囲】

【請求項1】

モータの振動現象の物理量を算出する算出方法であって、
前記モータが備えるステータコアの突極に作用する時間依存の加振力を、空間および時間の周波数ごとに変換して前記モータの振動現象の物理量を算出する算出工程を有し、
前記算出工程において、前記ステータコアの突極数に分割した領域のそれぞれに作用する値を前記加振力として用いる、ことを特徴とする算出方法。

【請求項2】

前記モータの振動現象の物理量は、円環次数ごとかつ時間周波数ごとに分解された振動値の分布、もしくは、前記振動値の分布の円環次数を前記ステータコアの突極数に分割した領域に対応して逆変換した振動値である、ことを特徴とする請求項１に記載の算出方法。

【請求項3】

前記モータの振動現象の物理量は、円環次数ごとかつ時間周波数ごとに分解された加振力の分布と、当該加振力の分布に対応して分解されたコンプライアンスの分布とをアダマール積算をして得られる振動変位の分布、もしくは、当該振動変位の分布の円環次数を前記ステータコアの突極数に分割した領域に対応して逆変換した時間依存の周波数の変位であることを特徴とする請求項１に記載の算出方法。

【請求項4】

前記モータの振動現象の物理量は、円環次数ごとかつ時間周波数ごとに分解された加振力の分布と、当該加振力の分布に対応して分解されたアクセラレンスの分布とをアダマール積算をして得られる振動加速度の分布、もしくは、当該振動加速度の分布の円環次数を前記ステータコアの突極数に分割した領域に対応して逆変換した時間依存の周波数の変位であることを特徴とする請求項１に記載の算出方法。

【請求項5】

前記算出工程において、前記ステータコアの突極数に分割した領域のそれぞれに作用する値のうち、径方向成分を使用することを特徴とする請求項１に記載の算出方法。

【請求項6】

前記算出工程において、前記加振力を径方向成分と接線方向成分に分離し、前記径方向成分と前記接線方向成分それぞれに基づいて振動現象の物理量を算出し、前記径方向成分と前記接線方向成分それぞれに基づいて算出された振動現象の物理量を、空間および時間の周波数ごとに合算することで、前記モータの振動現象の物理量を算出することを特徴とする請求項１に記載の算出方法。

【請求項7】

前記算出工程において、以下の式（１）が用いられる、ことを特徴とする請求項１に記載の算出方法。

【数1】

【請求項8】

請求項１～７のいずれか一項に記載の算出方法を実行する処理工程と、
前記処理工程にて算出された計算結果と、前記モータの速度および負荷の少なくとも一方に対応して予め規定された制御パラメータとに基づき、前記モータに対する印加電流を制御する制御工程と、
を有することを特徴とするモータ制御方法。

【請求項9】

モータの振動現象の物理量を算出する算出装置であって、
前記モータが備えるステータコアの突極に作用する時間依存の加振力を取得する取得手段と、
前記取得手段にて取得した加振力を用いて、空間および時間の周波数ごとに変換して前記モータの振動現象の物理量を算出する算出手段と
を有し、
前記取得手段は、前記ステータコアの突極数に分割した領域のそれぞれに作用する値を前記加振力として取得する、ことを特徴とする算出装置。

【請求項10】

コンピュータに、
モータが備えるステータコアの突極に作用する時間依存の加振力を、空間および時間の周波数ごとに変換して前記モータの振動現象の物理量を算出する算出工程を実行させ、
前記算出工程において、前記ステータコアの突極数に分割した領域のそれぞれに作用する値を前記加振力として用いる、ことを特徴とするプログラム。

【発明の詳細な説明】

【技術分野】

【0001】

本願発明は、モータの振動現象の物理量の算出方法、モータ制御方法、算出装置、およびプログラムに関する。

【背景技術】

【0002】

従来、様々な機器においてモータが利用されている。モータには様々な種類のものがあり、例えば、永久磁石同期電動機（ＰＭ（Ｐｅｒｍａｎｅｎｔ Ｍａｇｎｅｔ）モータ）が挙げられる。このようなモータの多くは、外周が円環状の積層コアで構成されている。モータ内の巻き線に電流を印加することで、電磁力を発生させ、その円周の接線成分をトルクとしてモータ出力を得ている。しかし、電磁力は接線方向の成分だけではなく、法線方向、すなわち、円環状の積層コアの径方向にも存在する。そのため、これらの成分が積層コアを径方向に変形させてしまい、振動や騒音などの有害な現象を引き起こしてしまう。

【0003】

上記のような現象を抑制するために、モータ内において発生する振動現象を特定することが求められている。公知の方法としては、例えば、有限要素法による構造解析や周波数応答解析などを用いて、振動現象のシミュレーションを行うことが知られている。有限要素法を用いる場合、物体を４面体や６面体などの単純な形状（要素）で分割し、その接点や辺に未知数を割り当てて剛性マトリクスを解くことなどが行われている。

【0004】

また、非特許文献１では、モータの空隙磁束密度を空間次数にて分解し、突極ごとを単位とする運動方程式モデルを立て、これを解くことで振動値（加速度）を算出する方法が提案されている。

【先行技術文献】

【非特許文献】

【0005】

【非特許文献1】鈴木照平ほか２名、「埋込磁石同期電動機（ＩＰＭＳＭ）における電磁共振条件の検討」，明電時報，Ｖｏｌ．３６７，Ｎｏ．２，ｐｐ．６－１１，２０２０年

【発明の概要】

【発明が解決しようとする課題】

【0006】

有限要素法を用いた構造解析は、複雑な形状を考慮できる一方、計算精度を確保するためには、要素分割を精細にする必要がある。そのため、未知数の個数を数百万規模で立てる行列を解く必要がある。よって、モータの構造に応じて計算量が増加し、その計算負荷により多くの所要時間を要するという問題がある。例えば、モータを用いた製品の開発期間や開発コストを抑制する観点からも、上記のような計算負荷は極力抑制することが求められている。

【0007】

また、非特許文献１の方法では、加振力を空隙磁束密度からマクスウェル応力法で求めている。モータの積層コアに作用する加振力は、区分的空隙に作用する応力の積分で求める必要がある。このとき、積層コアは有限の角度を持っていることから、径方向成分といえども突極周方向中央と端部では角度が異なる。そのため、積層コアを構成する突極が広いほど、つまり突極数が少ないほど誤差が大きくなるという問題がある。また、空隙磁束密度のデータ量は数百点あるため、計算量についても大きなものとなっていた。

【0008】

上記課題を鑑み、本願発明は、従来よりもさらに計算負荷を抑制し、計算の精度を向上させた、モータにおける振動現象の物理量の算出方法を提供することを目的とする。

【課題を解決するための手段】

【0009】

上記課題を解決するために本願発明は以下の構成を有する。すなわち、モータの振動現象の物理量を算出する算出方法は、
前記モータが備えるステータコアの突極に作用する時間依存の加振力を、空間および時間の周波数ごとに変換して前記モータの振動現象の物理量を算出する算出工程を有し、
前記算出工程において、前記ステータコアの突極数に分割した領域のそれぞれに作用する値を前記加振力として用いる。

【0010】

また、本願発明の別の形態は以下の構成を有する。すなわち、モータ制御方法は、
モータの振動現象の物理量を算出する算出方法を実行する処理工程と、
前記処理工程にて算出された計算結果と、前記モータの速度および負荷の少なくとも一方に対応して予め規定された制御パラメータとに基づき、前記モータに対する印加電流を制御する制御工程と、
を有し、
前記算出方法は、
前記モータが備えるステータコアの突極に作用する時間依存の加振力を、空間および時間の周波数ごとに変換して前記モータの振動現象の物理量を算出する算出工程を有し、
前記算出工程において、前記ステータコアの突極数に分割した領域のそれぞれに作用する値を前記加振力として用いる。

【0011】

また、本願発明の別の形態は以下の構成を有する。すなわち、モータの振動現象の物理量を算出する算出装置は、
前記モータが備えるステータコアの突極に作用する時間依存の加振力を取得する取得手段と、
前記取得手段にて取得した加振力を用いて、空間および時間の周波数ごとに変換して前記モータの振動現象の物理量を算出する算出手段と
を有し、
前記取得手段は、前記ステータコアの突極数に分割した領域のそれぞれに作用する値を前記加振力として取得する。

【0012】

また、本願発明の別の形態は以下の構成を有する。すなわち、プログラムは、
コンピュータに、
モータが備えるステータコアの突極に作用する時間依存の加振力を、空間および時間の周波数ごとに変換して前記モータの振動現象の物理量を算出する算出工程を実行させ、
前記算出工程において、前記ステータコアの突極数に分割した領域のそれぞれに作用する値を前記加振力として用いる。

【発明の効果】

【0013】

本願発明により、従来よりもさらに計算負荷を抑制し、計算の精度を向上させた、モータにおける振動現象を特定するための算出方法を提供することが可能となる。

【図面の簡単な説明】

【0014】

【図1】本願発明の一実施形態に係る算出方法を実行可能な装置構成の例を示す概略図。

【図2】本願発明の一実施形態に係る算出方法を適用可能なモータの概略構成を示す図。

【図3】本願発明の第１の実施形態に係る算出方法にて用いるコンプライアンスを示すグラフ図。

【図4】本願発明の第１の実施形態に係る算出方法にて用いるアクセラレンスを示すグラフ図。

【図5】本願発明の第１の実施形態に係る算出方法にて用いる磁界解析結果の例を示す図。

【図6】本願発明の第１の実施形態に係る算出方法における振幅値の例を示す図。

【図7】本願発明の第１の実施形態に係る算出方法における振幅値の例を示す図。

【図8】本願発明の第１の実施形態に係る算出方法における加速度の例を示す図。

【図9】本願発明の第１の実施形態に係る算出方法を用いた処理のフローチャート。

【図10】本願発明の第２の実施形態に係る算出方法にて用いるコンプライアンスを示すグラフ図。

【図11】本願発明の第２の実施形態に係る算出方法にて用いるアクセラレンスを示すグラフ図。

【図12】本願発明の第２の実施形態に係る算出方法にて用いる磁界解析結果の例を示す図。

【図13】本願発明の第２の実施形態に係る算出方法における加振力の径方向成分の例を示す図。

【図14】本願発明の第２の実施形態に係る算出方法における加振力の接線方向成分の例を示す図。

【図15】本願発明の第２の実施形態に係る算出方法における加速度の径方向成分の例を示す図。

【図16】本願発明の第２の実施形態に係る算出方法における加速度の接線方向成分の例を示す図。

【図17】本願発明の第２の実施形態に係る算出方法における加速度の例を示す図。

【図18】本願発明の第１の実施形態に係る算出方法における成分ごとの加速度の例を示す図。

【図19】本願発明の第２の実施形態に係る算出方法を用いた処理のフローチャート。

【発明を実施するための形態】

【0015】

以下、本願発明を実施するための形態について図面などを参照して説明する。なお、以下に説明する実施形態は、本願発明を説明するための一実施形態であり、本願発明を限定して解釈されることを意図するものではなく、また、各実施形態で説明されている全ての構成が本願発明の課題を解決するために必須の構成であるとは限らない。また、各図面において、同じ構成要素については、同じ参照番号を付すことにより対応関係を示す。

【0016】

＜第１の実施形態＞
以下、本願発明の第１の実施形態について説明を行う。本願発明に係る算出方法は、例えば、モータの振動現象のシミュレーションに用いることが可能であり、本実施形態では、シミュレーションにおいて利用可能な振動現象の物理量を算出する方法について説明する。

【0017】

［装置構成］
図１は、本実施形態に係る算出方法を実行可能な情報処理装置の全体構成の一例を示す概略構成図である。本実施形態に係る情報処理装置１は、例えば、ＰＣ（Ｐｅｒｓｏｎａｌ Ｃｏｍｐｕｔｅｒ）などが用いられてよく、その構成は特に限定するものではない。

【0018】

情報処理装置１は、ＣＰＵ（Ｃｅｎｔｒａｌ Ｐｒｏｃｅｓｓｉｎｇ Ｕｎｉｔ）１０、ＲＯＭ（Ｒｅａｄ Ｏｎｌｙ Ｍｅｍｏｒｙ）１１、ＲＡＭ（Ｒａｎｄｏｍ Ａｃｃｅｓｓ Ｍｅｍｏｒｙ）１２、ＨＤＤ（Ｈａｒｄ Ｄｉｓｋ Ｄｒｉｖｅ）１３、入力装置１４、表示装置１５、および通信装置１６を含んで構成される。ＣＰＵ１０は、情報処理装置１全体の制御を司る部位であり、例えば、ＨＤＤ１３に格納されたプログラムを読み出して実行することで各種機能を実現してよい。ＲＯＭ１１は、不揮発性の記憶領域である。ＲＡＭ１２は、揮発性の記憶領域であり、一時的なデータの保存場所として用いられる。ＨＤＤ１３は、不揮発性の記憶領域であり、各種プログラムやデータが記憶、管理される。入力装置１４は、外部からの入力を受け付ける部位であり、例えば、マウスやキーボードなどから構成される。表示装置１５は、各種情報を表示するための部位であり、例えば、液晶ディスプレイなどが該当する。なお、入力装置１４と表示装置１５が一体となったタッチパネルディスプレイが用いられてもよい。通信装置１６は、外部装置（不図示）とネットワーク（不図示）を介して通信するための部位である。ここでの通信は、有線／無線は問わず、また、通信規格なども特に限定するものではない。

【0019】

情報処理装置１は、後述する本実施形態に係る算出方法を実行するためのプログラムの他、電磁界解析を行うための汎用の電磁界ソフトウェアを実行可能な構成であってよい。更には、本実施形態に係る算出方法で用いられる、変換係数（コンプライアンスやアクセラレンス）を算出するためのプログラムを実行可能であってよい。なお、これらのソフトウェアやプログラムは、必ずしも後述する算出方法を実行する装置と同じ装置で実行される必要はない。例えば、上記のようなプログラムが別個の装置で実行された後、その処理結果が本実施形態に係る算出方法が実行される装置に入力されるような構成であってもよい。

【0020】

［モータ概略］
図２は、本実施形態に係るモータ２００の概略構成を示す図である。本実施形態に係る算出方法を適用可能なモータ２００として、永久磁石同期電動機（ＰＭモータ）が挙げられる。なお、ここでは説明を省略するが、モータ２００には、モータ２００の動作を制御する制御装置（不図示）や、モータ２００により出力されるトルクを伝達するための伝達部（不図示）などが接続されてよい。また、ＰＭモータを動作させるためのインバータ（不図示）や、モータ２００の回転子の位置を検出するための位置センサ（不図示）などが備えられてよい。

【0021】

図２（ａ）は、モータ２００の断面の概略構成を示す。モータ２００は、固定子（ステータ）２０１と、回転子（ロータ）２０２とを含んで構成される。ロータ２０２は、その円周に沿って、複数の永久磁石２０３が設けられる。ここでは、モータ２００として、永久磁石２０３がロータ２０２の内部に設けられたＩＰＭ（Ｉｎｔｅｒｉｏｒ Ｐｅｒｍａｎｅｎｔ Ｍａｇｎｅｔ：磁石埋め込み式）モータを例に挙げ、６極の構成例を示している。しかし、この構成に限定するものではなく、ロータの表面に磁石を張り合わせたＳＰＭ（Ｓｕｒｆａｃｅ Ｐｅｒｍａｎｅｎｔ Ｍａｇｎｅｔ）モータが用いられてもよい。ステータ２０１とロータ２０２の間には、エアギャップが存在する。また、図２（ａ）では不図示であるが、ステータ２０１には、モータの駆動に係る電磁力を発生させるための巻き線等によるコイルが設けられる。

【0022】

図２（ｂ）は、図２（ａ）に示すモータ２００の構成のうち、ステータ２０１の固定子コア（ステータコア）にのみ着目した概略図である。ステータ２０１のステータコアは、複数のティース２０４を含んで構成され、ここでは、１２個のティース２０４ａ～２０４ｌから構成される例を示している。以下の説明において、ティースを包括的に説明する場合には添え字（ａ～ｌ）を省略して示し、個別に説明を要する場合には、添え字を付して説明する。本実施形態では、複数のティース２０４について、加振力に応じた領域分けを行う。図２（ｂ）の場合、１２個のティース２０４それぞれに対応したロータ周りの周方向の領域を対象として、１２の領域に分割した例を示す。

【0023】

［算出方法］
本実施形態では、図２で示したような環状物であるステータ２０１に対し、バネマス系運動方程式として、振動現象を示す物理量の算出を行う。ステータ２０１は、振動モード次数（以下、円環次数とも称する）ごとにバネ定数（固有値）が異なる。そのため、周波数ごと、円環次数ごとに方程式を立てる。円環次数としては、図２に示したステータ２０１の構成の場合、０次、２次、４次、６次が挙げられる。

【0024】

ここで、従来のバネマス系運動方程式として、以下の式（１）～式（３）が知られている。

【0025】

【数1】

【0026】

ｕ：突極数
ｉ：円環次数（振動モード次数）
ｍ：ステータの全重量
ｒ：ステータ外周のラジアル変位
ω_０ｉ：ｉ次（円環次数）の固有値角速度（周波数の２π倍）
ω’：加振力の角速度
Ｂ：加振力密度の径方向成分振幅（添え字ｃはｃｏｓ成分、添え字ｓはｓｉｎ成分）
ｊ：虚数単位
ｇ：加振力の時間ｔの関数

【0027】

非特許文献１では、空隙磁束密度を用いて運動方程式を立てている。このとき、ステータの突極数が少ない場合などには突極間の開口部（ロータ周りにおける突極間の隙間）が生じ、その結果、加振力は各突極に集約されることが想定される。加振力が集約することにより、空隙磁束密度の分布には誤差が生じ、算出結果の精度が低下する恐れがある。これは、ステータの突極の数が少なくなるに従って、その誤差が顕著になり得る。

【0028】

そこで、本実施形態では、図２（ｂ）に示したような分割の範囲における加振力を接点力法で求める。つまり、加振力がステータ２０１のステータコアに作用する加振力分布に相似であると仮定し、上記の式（１）を以下の式（４）のように修正する。

【0029】

【数2】

【0030】

δ：円環次数ごとの補正係数
Ｆ：加振力の径方向成分振幅（添え字ｃはｃｏｓ成分、添え字ｓはｓｉｎ成分）
Ｃｍｐ：コンプライアンス（＝変位／加振力）
Ａｃｌ：アクセラレンス（＝加速度／加振力）

【0031】

つまり、式（１）では、加振力密度の径方向成分振幅Ｂを用いているのに対し、式（４）では、加振力の径方向成分振幅Ｆを用いている。上記の式（４）のように修正することで、複数の突極における加振力の集約を考慮し、突極の形状に依存せずに計算が可能となる。

【0032】

上記の式（４）に関し、周期振動となる十分な時間が経過した特解は以下の式のようになる。

【0033】

【数3】

【0034】

【数4】

【0035】

【数5】

【0036】

【数6】

【0037】

【数7】

【0038】

変換係数であるコンプライアンスＣｍｐおよびアクセラレンスＡｃｌは、構造解析の結果、および、上記の式（８）および式（９）を用いて導出する。なお、変換係数の計算数はごく少数であり、予め算出しておくことで、以降の計算にて再利用可能である。そのため、各変換係数は、予め導出され、情報処理装置１の記憶部等に保持されていてよい。

【0039】

円環次数ごとの固有値である角速度ω_０ｉは、構造解析（固有値解析）の結果に基づいて導出する。また、疑似質量δ_ｉｍは、角速度ω’（周波数ｆ’＝ω’／２π）がごく小さい状態における構造解析（周波数応答解析）の結果に基づいて、以下の式（１０）、式（１１）を用いて導出する。

【0040】

【数8】

【0041】

そして、式（１０）、式（１１）で求めた値を、式（８）、式（９）に代入することで、変換係数であるコンプライアンスＣｍｐおよびアクセラレンスＡｃｌを導出する。なお、上記の式において、ω’とω_０とが近似した際に、コンプライアンスＣｍｐやアクセラレンスＡｃｌに特異点となる値が生じ得る。その場合には、上記式に、モーダル減衰比を適用して調整することで算出を行ってもよい。モーダル減衰比の具体的な値については、振動解析ソフトなどで用いられている一般的な値を用いてよく、特に限定するものではない。また、特異点を回避するためにそのほかの公知のパラメータを含めてもよい。

【0042】

図３は、上記式により導出された、円環次数ごとの周波数とコンプライアンスの関係を示すグラフである。図３において、横軸は周波数［Ｈｚ］を示し、縦軸はコンプライアンス［任意単位］を示す。ここでは、円環次数として、０次、２次、４次、６次の値をそれぞれ示している。

【0043】

図４は、上記式により導出された、円環次数ごとの周波数とアクセラレンスの関係を示すグラフである。図４において、横軸は周波数［Ｈｚ］を示し、縦軸はアクセラレンス［任意単位］を示す。ここでは、円環次数として、０次、２次、４次、６次の値をそれぞれ示している。

【0044】

［処理フロー］
図９は、本実施形態に係るモータ２００の振動現象に係る物理量を算出するための算出処理のフローチャートである。本処理は、情報処理装置１により実行され、例えば、情報処理装置１が備えるＣＰＵ１０が本実施形態に係る処理を実現するためのプログラムをＨＤＤ１３から読み出して実行することにより実現される。

【0045】

Ｓ９０１にて、情報処理装置１は、円環次数ごとの周波数に対するコンプライアンス（変換係数Ｃｍｐ）およびアクセラレンス（変換係数Ａｃｌ）を入力する。上述したように、コンプライアンスおよびアクセラレンスは、構造解析の結果に基づき、式（８）～式（１１）を用いて導出されて保持されておき、情報処理装置１が本処理フローの際に参照するように構成されていてよい。

【0046】

Ｓ９０２にて、情報処理装置１は、モータ２００のステータ２０１（ステータコア）が備える突極ごとの時間に対する加振力を入力する。加振力は、電磁界解析結果から取得されてよい。電磁界解析結果は、汎用の磁界解析ソフトなどを用いて２次元の電磁界解析を行うことで得られてよい。図５は、電磁界解析結果から得られる加振力の例を示す。図５において、横軸は電気角［ｄｅｇ］を示し、縦軸は加振力［任意単位］を示す。ここでは、ステータ２０１が１８個の突極（ティース）を備えるモータにおける加振力の解析結果の例を示す。以下の式（１２）は、突極ごとの時間に対する加振力（実数）の行列式を示す。つまり、式（１２）は、突極と時間に応じて分解した加振力の分布を示す行列式である。

【0047】

【数9】

【0048】

ｆ：加振力
ｔ：時間（ｔ_０～ｔ_ｍ－１）
Ｔ：突極番号（Ｔ_１～Ｔ_ｕ）
ｕ：突極数
ｍ：時間分割数

【0049】

Ｓ９０３にて、情報処理装置１は、Ｓ９０２にて入力された突極ごとの時間に対する加振力に対して２次元ＦＦＴ（Ｆａｓｔ Ｆｏｕｒｉｅｒ Ｔｒａｎｓｆｏｒｍ）を適用することで、円環次数ごとの周波数に対する加振力へ変換する。ＦＦＴについては公知の方法を用いてよく、ここでの詳細な説明は省略する。以下の式（１３）は、式（１２）に対して２次元ＦＦＴ処理を行う適用した後の、円環次数ごとの周波数に対する加振力（複素数）の行列式を示す。つまり、式（１３）は、円環次数と周波数に応じて分解した加振力の分布を示す行列式である。

【0050】

【数10】

【0051】

ｃ：円環次数（ｃ_－ｕｈ～ｃ_ｕｈ；ｕｈ＝ｕ／２；正負の符号は回転方向を示す）
ｚ：周波数（ｚ_０～ｚ_ｍｈ；ｍｈ＝ｍ／２）

【0052】

図６は、円環次数ごとの周波数に対する振幅値を示すグラフの例を示す。図６において、横軸は周波数［Ｈｚ］を示し、縦軸は振幅値［任意単位］を示す。

【0053】

Ｓ９０４にて、情報処理装置１は、加振力とコンプライアンスＣｍｐの乗算による円環次数ごとの周波数に対する変位（振動変位）を算出する。コンプライアンスＣｍｐは、図３に示したような値が用いられる。以下の式（１４）は、コンプライアンスＣｍｐを示す行列式である。また、以下の式（１５）は、式（１３）と式（１４）から変位を算出するための式を示す。なお、式（１５）において、算術記号〇は、アダマール積を示す。つまり、式（１３）に示す行列式と、式（１４）に示す行列式をアダマール積算することで、円環次数と周波数に応じて分解された振動値（物理量）を求める。

【0054】

【数11】

【0055】

【数12】

【0056】

図７は、Ｓ９０４の処理結果として得られる円環次数ごとの周波数に対する振幅値を示すグラフの例を示す。図７において、横軸は周波数［Ｈｚ］を示し、縦軸は振幅値［任意単位］を示す。

【0057】

Ｓ９０５にて、情報処理装置１は、加振力とアクセラレンスＡｃｌの乗算による円環次数ごとの周波数に対する加速度（振動加速度）を算出する。アクセラレンスＡｃｌは、図４に示したような値が用いられる。ここでの算出式は、Ｓ９０４にて用いたコンプライアンスの場合の式（１５）と同様である。ここでは、図４に示すような値に合わせて、式（１４）におけるパラメータが異なる。算出結果を示すグラフについても、図７と同様の構成となる。

【0058】

Ｓ９０６にて、情報処理装置１は、Ｓ９０４およびＳ９０５にて算出された円環次数ごとの周波数に対する振動値（変位および加速度）に対し、実数２次元逆ＦＦＴを適用することで、突極ごとの時間に対する変位および加速度に変換する。実数２次元逆ＦＦＴについては、Ｓ９０３の２次元ＦＦＴに対応して公知の方法を用いてよく、ここでの詳細な説明は省略する。以下の式（１６）は、実数２次元逆ＦＦＴ後の振動値を示す行列式である。つまり、式（１６）は、突極と時間に応じて分解した振動値の分布を示す行列式である。

【0059】

【数13】

【0060】

ｘ：振動値（加速度、または、変位）

【0061】

Ｓ９０７にて、情報処理装置１は、Ｓ９０６にて得られた突極ごとの時間に対する変位および加速度それぞれに対し、実数１次元ＦＦＴを適用することで、突極ごとの周波数に対する変位および加速度に変換する。実数１次元ＦＦＴについては、公知の方法を用いてよく、ここでの詳細な説明は省略する。以下の式（１７）は、実数１次元ＦＦＴ後の振動値を示す行列式である。つまり、式（１７）は、突極と周波数に応じて分解した振動値の分布を示す行列式である。

【0062】

【数14】

【0063】

Ｓ９０８にて、情報処理装置１は、Ｓ９０７で得られた振幅値を用いて、突極ごとの周波数に対する変位および加速度それぞれの最大値を算出する。以下の式（１８）は、Ｓ９０７にて得られた値に基づいて、周波数ごとの最大振幅、すなわち、変位の最大値および加速度の最大値を求めるための式を示す。なお、式（１８）におけるｍａｘ関数は、複数の引数のうちの最大値を求める関数である。

【0064】

【数15】

【0065】

図８は、Ｓ９０８の処理結果として得られる、評価値（ここでは、加速度）と周波数との関係を示すグラフの例を示す。図８において、横軸は周波数［Ｈｚ］を示し、縦軸は評価値［任意単位］を示す。図８では本手法に対する比較対象として、有限要素法を用いた構造解析により得られた結果を示している。図８に示すように、本手法による評価値は、構造解析とほぼ同等の値を示している。

【0066】

次に、本実施形態に係る手法と、従来の手法による計算時間との差異について説明する。ここで、従来におけるモータの振動解析を行う手法の例として、有限要素法を用いた構造解析を使用する手法（以下、「従来手法１」とも称する）と、空隙磁束密度基準による手法（非特許文献１）（以下、「従来手法２」とも称する）を挙げて説明する。

【0067】

モータの振動解析を行う場合、従来手法１では、磁界解析と構造解析の２つの工程を要する。また、従来手法２および本手法では、磁界解析と振動計算の２つの工程を要する。いずれの手法においても、磁界解析に要する時間は同等であるため、ここではもう一方の工程に着目して説明する。なお、解析の事前準備に要する時間についても省略して説明する。なお、以下に示す数値は一例であり、解析環境や解析対象により変動が生じ得るが、ここではすべて同じ解析環境下で、同じ解析対象を用いた例を示す。

【0068】

従来手法１の構造解析については、１回の解析当たり、４８００～６０００秒程度の計算時間を要する。従来手法２の振動計算では、１回の解析当たり、５０秒程度の計算時間を要する。本手法の振動計算では、１回の解析当たり、２秒程度の計算時間を要する。なお、いずれの手法でも必要となる磁界解析は、１５０秒程度を要する。

【0069】

図８にて示したように、従来手法１と本手法とは、ほぼ同等の精度を得られる。その一方、本手法の振動計算に要する時間は、従来手法１の構造解析に要する時間と比較して、極めて小さい値となる。上記の例では、計算時間に数千倍の差がある。

【0070】

また、本手法の振動計算に要する時間は、従来手法２の振動計算に要する時間と比較しても、小さい値となる。上記の例では、計算時間に数十倍の差がある。

【0071】

モータの振動解析を行う前提として、様々な条件下（例えば、モータ磁気回路の形状や電流波形、回転数、トルクなどに係る条件）での解析を要するため、解析においては複数回の計算（例えば、数千～数万）を複数回行う必要がある。例えば、モータの低振動、低騒音を実現するために、条件を切り替えてより多くの計算を行って最適化を行う必要がある。このようなモータを用いた製品の開発における状況下では、上記のような計算時間の差は非常に大きく、その有用性はより多くの振動解析を行う必要のある製品開発に対してより顕著となる。

【0072】

更に、従来手法２と本手法との解析精度の観点から説明する。非特許文献１における課題について上述したように、従来手法２では、モータの積層コアに作用する加振力は、区分的空隙に作用する応力の積分で求める必要がある。このとき、モータの積層コアを行使する突極が広いほど、つまり突極数が少ないほど誤差が大きくなるという問題がある。そのため、適用するモータの構成によっては、従来手法２では、従来手法１のような精度を得ることができない。つまり、従来手法２では、従来手法１と同等の精度を得ることができる本手法ほどの精度を得ることができない。また、従来手法２では、振動を広い周波数帯域で見る必要があり、これは、数十次程度の高次周波数成分が影響してしまう。つまり、低振動の成分については、適切に測定することができない。そのため、例えば、モータの開発において低振動や低騒音の実現のために振動解析を行う場合に、上記のような最適化（条件調整）を行うための解析には、従来手法２を用いることが困難となる。

【0073】

以上の処理により、本実施形態では、モータの振動現象を、例えば、有限要素法を用いた手法などの従来手法よりも少ない計算量により、同等の精度の算出が可能となる。また、本実施形態に係る手法により、モータにおける加振力の突極における集約を考慮して、振動値を算出することができる。更には、本実施形態に係る手法では、モータの突極の形状に関わらず、例えば、空隙磁束密度基準による手法などの従来手法よりも精度よく振動値を算出することができる。なお、図９に示す処理フローの並びはこれに限定するものではなく、一部の処理が入れ替わってもよい。例えば、Ｓ９０４の変位の算出と、Ｓ９０５の加速度の算出は逆であってよい。また、変位と加速度のいずれか一方のみを算出したい場合には、一方の算出工程が省略されてもよい。

【0074】

本実施形態に係るモータの振動現象の解析は、例えば、自動車の開発、より具体的にはＥＰＳ（Ｅｌｅｃｔｒｉｃ Ｐｏｗｅｒ Ｓｔｅｅｒｉｎｇ）装置の開発の際に適用可能である。上述したように、本願発明に係る手法では、従来手法よりも少ない計算負荷で、従来と同等の計算精度を得ることができる。また、例えば、低振動や低騒音のモータ開発における条件最適化のための探索などにも利用することができる。そして、本手法により、装置の設計期間の短縮や工程削減によるコストの削減、試験の簡易化などを実現することが可能となる。更には、上記効果を踏まえ、製品自体の品質向上なども促進することが可能となる。

【0075】

＜第２の実施形態＞
本願発明の第２の実施形態について説明する。なお、第１の実施形態と重複する構成については説明を省略し、差分に着目して説明する。

【0076】

第１の実施形態では、磁界解析にて求められた加振力と、構造解析にて求められた変換係数を用いて、振動現象を示す物理量を算出していた。本実施形態では、加振力を更にステータコアの径方向および接線方向の各成分に分解して用いる。そして、径方向成分と、接線方向成分に対して、第１の実施形態で説明した手法と同様に予め変換係数（コンプライアンスおよびイナータンス）を算出しておき、それらの値を用いて振動現象を表す物理量を算出する。

【0077】

［算出方法］
本実施形態では、第１の実施形態と同様、図２で示したような環状物であるステータ２０１に対し、バネマス系運動方程式として、振動現象を示す物理量の算出を行う。バネマス系運動方程式は、上記の式（１）～式（３）を用いる。なお、以下では解析単位をブロックとも称する。

【0078】

次に、本実施形態では、図２（ｂ）に示したような分割の範囲における加振力の径方向成分および接線方向成分それぞれを接点力法で求める。

【0079】

（径方向成分）
まずは、径方向成分について示す。加振力がステータ２０１のステータコアに作用する加振力分布に相似であると仮定し、上記の式（１）を以下の式（１９）のように修正する。

【0080】

【数16】

【0081】

δ_Ｒ：径方向成分における円環次数ごとの補正係数（ブロック標本化比例定数）
Ｆ_Ｒ：加振力の径方向成分振幅（添え字ｃはｃｏｓ成分、添え字ｓはｓｉｎ成分）
ｔ：時間
ｅ：自然対数の底
ｍ：ステータ全質量
ｒ：径方向変位
ω_０：径方向固有値角速度
ω’：角速度
ζ：ダンピング係数
Ｃｍｐ：コンプライアンス
Ｍｏｖ：モビリティ
Ｉｎｒ：イナータンス

【0082】

つまり、式（１）では、加振力密度の径方向成分振幅Ｂを用いているのに対し、式（１９）では、加振力の径方向成分振幅Ｆ_Ｒを用いている。上記の式（１９）のように修正することで、複数の突極における径方向成分の加振力の集約を考慮し、突極の形状に依存せずに計算が可能となる。

【0083】

上記の式（１９）に関し、周期振動となる十分な時間が経過した特解は以下の式のようになる。

【0084】

【数17】

【0085】

【数18】

【0086】

【数19】

【0087】

【数20】

【0088】

【数21】

【0089】

【数22】

【0090】

【数23】

【0091】

【数24】

【0092】

【数25】

【0093】

変換係数である径方向成分のコンプライアンスＣｍｐ_Ｒｉ、径方向成分のモビリティＭｏｖ_Ｒｉ、および径方向成分のイナータンスＩｎｒ_Ｒｉは、構造解析の結果、および、上記の式（２６）～式（２８）を用いて導出する。イナータンスは、第１の実施形態におけるアクセラレンスに対応する。なお、変換係数の計算数はごく少数であり、予め算出しておくことで、以降の計算にて再利用可能である。そのため、各変換係数は、予め導出され、情報処理装置１の記憶部等に保持されていてよい。

【0094】

円環次数ごとの固有値である径方向成分の角速度ω_Ｒ０ｉは、構造解析（固有値解析）の結果に基づいて導出する。また、径方向成分の疑似質量δ_Ｒｉｍは、角速度ω’（周波数ｆ’＝ω’／２π）がごく小さい状態における構造解析（周波数応答解析）の結果に基づいて、以下の式（２９）、式（３０）を用いて導出する。

【0095】

【数26】

【0096】

【数27】

【0097】

そして、式（２９）、式（３０）で求めた値を、式（２６）～式（２８）に代入することで、変換係数である径方向成分のコンプライアンスＣｍｐ_Ｒｉ、径方向成分のモビリティＭｏｖ_Ｒｉ、および径方向成分のイナータンスＩｎｒ_Ｒｉを導出する。

【0098】

（接線方向成分）
次に接線方向成分について示す。接線方向についても、径方向成分と同様に、加振力がステータ２０１のステータコアに作用する加振力分布に相似であると仮定し、上記の式（１）を以下の式（３１）のように修正する。

【0099】

【数28】

【0100】

δ_Ｔ：接線方向成分における円環次数ごとの補正係数（ブロック標本化比例定数）
Ｆ_Ｔ：加振力の接線方向成分振幅（添え字ｃはｃｏｓ成分、添え字ｓはｓｉｎ成分）

【0101】

つまり、式（３１）では、加振力の接線方向成分振幅Ｆ_Ｔを用いている。上記の式（３１）のように修正することで、複数の突極における接線方向の加振力の集約を考慮し、突極の形状に依存せずに計算が可能となる。

【0102】

上記の式（３１）に関し、周期振動となる十分な時間が経過した特解は以下の式のようになる。

【0103】

【数29】

【0104】

【数30】

【0105】

【数31】

【0106】

【数32】

【0107】

【数33】

【0108】

【数34】

【0109】

【数35】

【0110】

【数36】

【0111】

【数37】

【0112】

変換係数である接線方向成分のコンプライアンスＣｍｐ_Ｔｉ、接線方向成分のモビリティＭｏｖ_Ｔｉ、および接線方向成分のイナータンスＩｎｒ_Ｔｉは、構造解析の結果、および、上記の式（３８）～式（４０）を用いて導出する。イナータンスは、第１の実施形態におけるアクセラレンスに対応する。なお、変換係数の計算数はごく少数であり、予め算出しておくことで、以降の計算にて再利用可能である。そのため、各変換係数は、予め導出され、情報処理装置１の記憶部等に保持されていてよい。

【0113】

円環次数ごとの固有値である接線方向成分の角速度ω_Ｔ０ｉは、構造解析（固有値解析）の結果に基づいて導出する。また、接線方向成分の疑似質量δ_Ｔｉｍは、角速度ω’（周波数ｆ’＝ω’／２π）がごく小さい状態における構造解析（周波数応答解析）の結果に基づいて、以下の式（４１）、式（４２）を用いて導出する。

【0114】

【数38】

【0115】

【数39】

【0116】

そして、式（４１）、式（４２）で求めた値を、式（３８）～式（４０）に代入することで、変換係数である接線方向成分のコンプライアンスＣｍｐ_Ｔｉ、接線方向成分のモビリティＭｏｖ_Ｔｉ、および接線方向成分のイナータンスＩｎｒ_Ｔｉを導出する。

【0117】

図１０は、上記式により導出された、径方向および接線方向それぞれの円環次数ごとの周波数とコンプライアンスの関係を示すグラフである。図１０において、横軸は周波数［Ｈｚ］を示し、縦軸はコンプライアンス［任意単位］を示す。ここでは、円環次数の例として、０次、２次、６次の値をそれぞれ示している。

【0118】

図１１は、上記式により導出された、径方向および接線方向それぞれの円環次数ごとの周波数とアクセラレンスの関係を示すグラフである。図１１において、横軸は周波数［Ｈｚ］を示し、縦軸はアクセラレンス［任意単位］を示す。ここでは、円環次数の例として、０次、２次、６次の値をそれぞれ示している。

【0119】

［処理フロー］
図１９は、本実施形態に係るモータ２００の振動現象に係る物理量を算出するための算出処理のフローチャートである。本処理は、情報処理装置１により実行され、例えば、情報処理装置１が備えるＣＰＵ１０が本実施形態に係る処理を実現するためのプログラムをＨＤＤ１３から読み出して実行することにより実現される。

【0120】

Ｓ１９０１にて、情報処理装置１は、径方向および接線方向それぞれに対して、円環次数ごとの周波数に対するコンプライアンス（変換係数Ｃｍｐ）およびイナータンス（変換係数Ｉｎｒ）を入力する。上述したように、コンプライアンスおよびイナータンスは、構造解析の結果に基づき、式（２６）、式（２８）、式（３８）、式（４０）を用いて導出されて保持されておき、情報処理装置１が本処理フローの際に参照するように構成されていてよい。径方向の変換係数を式（４３）に示し、接線方向の変換係数を式（４４）に示す。

【0121】

【数40】

【0122】

【数41】

【0123】

Ｙ：係数
Ｒ＿ｃ：径方向の円環次数（ｃ_－ｕｈ～ｃ_ｕｈ；ｕｈ＝ｕ／２；正負の符号は回転方向を示す）
Ｔ＿ｃ：接線方向の円環次数（ｃ_－ｕｈ～ｃ_ｕｈ；ｕｈ＝ｕ／２；正負の符号は回転方向を示す）
ｚ：周波数（ｚ_０～ｚ_ｍｈ；ｍｈ＝ｍ／２）

【0124】

Ｓ１９０２にて、情報処理装置１は、モータ２００のステータ２０１（ステータコア）が備える突極ごとの時間に対する加振力を入力する。加振力は、電磁界解析結果から取得されてよい。電磁界解析結果は、汎用の磁界解析ソフトなどを用いて２次元の電磁界解析を行うことで得られてよい。図１２は、電磁界解析結果から得られる径方向成分および接線方向成分それぞれの加振力の例を示す。図１２において、横軸は電気角［ｄｅｇ］を示し、縦軸は加振力［任意単位］を示す。ここでは、ステータ２０１が１８個の突極（ティース）、すなわちブロックを備えるモータにおける、径方向成分および接線方向成分の加振力の解析結果の例を示す。本実施形態では、電磁界解析により、モータ２００における加振力として、径方向成分の加振力と、接線方向成分の加振力をそれぞれ取得して用いる。以下の式（４５）は、突極ごとの時間に対する径方向の加振力（実数）の行列式を示す。また、以下の式（４６）は、突極ごとの時間に対する接線方向の加振力（実数）の行列式を示す。つまり、式（４５）、式（４６）は、突極と時間に応じて分解した、径方向成分および接線方向成分それぞれの加振力の分布を示す行列式である。

【0125】

【数42】

【0126】

【数43】

【0127】

ｆ_Ｒ：加振力（径方向成分）
ｆ_Ｔ：加振力（接線方向成分）
ｔ：時間（ｔ_０～ｔ_ｍ－１）
Ｔ：突極番号（Ｔ_１～Ｔ_ｕ）
ｕ：突極数
ｍ：時間分割数

【0128】

Ｓ１９０３にて、情報処理装置１は、Ｓ１９０２にて入力された突極ごとの時間に対する加振力に対して２次元ＦＦＴ（Ｆａｓｔ Ｆｏｕｒｉｅｒ Ｔｒａｎｓｆｏｒｍ）を適用することで、円環次数ごとの周波数に対する加振力へ変換する。ＦＦＴについては公知の方法を用いてよく、ここでの詳細な説明は省略する。以下の式（４７）は、式（４５）および式（４６）に対して、２次元ＦＦＴ処理を行う適用した後の、径方向成分および接線方向における円環次数ごとの周波数に対する加振力（複素数）の行列式を示す。つまり、式（４７）は、円環次数と周波数に応じて分解した加振力の分布を示す行列式である。

【0129】

【数44】

【0130】

図１３は、円環次数ごとの周波数に対する加振力の径方向成分を示すグラフの例を示す。図１３において、横軸は周波数［Ｈｚ］を示し、縦軸は加振力［任意単位］を示す。

【0131】

図１４は、円環次数ごとの周波数に対する加振力の接線方向成分を示すグラフの例を示す。図１４において、横軸は周波数［Ｈｚ］を示し、縦軸は加振力［任意単位］を示す。

【0132】

Ｓ１９０４にて、情報処理装置１は、Ｓ１９０３にて得られた径方向および接線方向の加振力と、Ｓ１９０１にて入力された径方向および接線方向それぞれのコンプライアンスＣｍｐの乗算による円環次数ごとの周波数に対する変位（振動変位）を算出する。コンプライアンスＣｍｐ（Ｃｍｐ_Ｒｉ、Ｃｍｐ_Ｔｉ）は、上述したように、径方向成分と接線方向成分それぞれに対応して設けられる。以下の式（４８）は、式（４３）、式（４４）、式（４７）から変位を算出するための式を示す。なお、式（４８）において、算術記号〇は、アダマール積を示す。つまり、式（４７）に示す行列式と、式（４３）または式（４４）に示す行列式とをアダマール積算することで、径方向成分および接線方向成分それぞれに対し、円環次数と周波数に応じて分解された振動値（物理量）を求める。

【0133】

【数45】

【0134】

Ｓ１９０５にて、情報処理装置１は、Ｓ１９０３にて得られた径方向および接線方向の加振力と、Ｓ１９０１にて入力された径方向および接線方向それぞれのイナータンスＩｎｒの乗算による円環次数ごとの周波数に対する加速度（振動加速度）を算出する。ここでの算出式は、Ｓ１９０４にて用いたコンプライアンスの場合の式（４８）と同様である。

【0135】

図１５は、Ｓ１９０５の処理結果として得られる円環次数ごとの周波数に対する加速度の径方向成分を示すグラフの例を示す。図１５において、横軸は周波数［Ｈｚ］を示し、縦軸は加速度［任意単位］を示す。

【0136】

図１６は、Ｓ１９０５の処理結果として得られる円環次数ごとの周波数に対する加速度の接線方向成分を示すグラフの例を示す。図１６において、横軸は周波数［Ｈｚ］を示し、縦軸は加速度［任意単位］を示す。

【0137】

Ｓ１９０６にて、情報処理装置１は、Ｓ１９０４にて算出した径方向成分に起因する変位と、接線方向成分に起因する変位とを、空間および時間の周波数ごとに合算することで、両方の成分を考慮した総合的な変位を求める。同様に、情報処理装置１は、Ｓ１９０５にて算出した径方向成分に起因する加速度と、接線方向成分に起因する加速度とを、空間および時間の周波数ごとに合算することで、総合的な加速度を求める。以下の式（４９）は、式（４８）から総合的な振動値（変位または加速度）を算出するための式を示す。

【0138】

【数46】

【0139】

ｘ_Ｒ：振動値（径方向成分由来の変位または加速度）
ｘ_Ｔ：振動値（接線方向成分由来の変位または加速度）
ｘ：振動値（変位または加速度）

【0140】

Ｓ１９０７にて、情報処理装置１は、Ｓ１９０６にて算出された円環次数ごとの周波数に対する振動値（変位および加速度）に対し、実数２次元逆ＦＦＴを適用することで、突極ごとの時間に対する変位および加速度に変換する。実数２次元逆ＦＦＴについては、Ｓ１９０３の２次元ＦＦＴに対応して公知の方法を用いてよく、ここでの詳細な説明は省略する。以下の式（５０）は、実数２次元逆ＦＦＴ後の振動値を示す行列式である。つまり、式（５０）は、突極と時間に応じて分解した振動値の分布を示す行列式であり、第１の実施形態にて示した式（１６）と同様である。

【0141】

【数47】

【0142】

Ｓ１９０８にて、情報処理装置１は、Ｓ１９０７にて得られた突極ごとの時間に対する変位および加速度それぞれに対し、実数１次元ＦＦＴを適用することで、突極ごとの周波数に対する変位および加速度に変換する。実数１次元ＦＦＴについては、公知の方法を用いてよく、ここでの詳細な説明は省略する。以下の式（５１）は、実数１次元ＦＦＴ後の振動値を示す行列式である。つまり、式（５１）は、突極と周波数に応じて分解した振動値の分布を示す行列式であり、第１の実施形態にて示した式（１７）と同様である。

【0143】

【数48】

【0144】

Ｓ１９０９にて、情報処理装置１は、Ｓ１９０８で得られた振幅値を用いて、突極ごとの周波数に対する変位および加速度それぞれの最大値を算出する。以下の式（５１）は、Ｓ１９０８にて得られた値に基づいて、周波数ごとの最大振幅、すなわち、変位の最大値および加速度の最大値を求めるための式を示し、第１の実施形態にて示した式（１８）と同様である。

【0145】

【数49】

【0146】

図１７は、Ｓ１９０９の処理結果として得られる、評価値（ここでは、加速度）と周波数との関係を示すグラフの例を示す。図１７において、横軸は周波数［Ｈｚ］を示し、縦軸は評価値［任意単位］を示す。図１７では本手法に対する比較対象として、有限要素法を用いた構造解析により得られた結果を示している。図１７に示すように、本手法による評価値は、構造解析とほぼ同等の値を示している。

【0147】

図１８は、図１７にて示している本手法による評価値を、径方向成分と接線方向成分それぞれに分けて示したグラフの例を示す。図１８において、横軸は周波数［Ｈｚ］を示し、縦軸は評価値［任意単位］を示す。図１７および図１８によると、径方向成分と接線方向成分の両方を考慮して算出を行うことで、モータの振動現象の物理量をより精度良く導くことができている。

【0148】

以上、本実施形態の構成により、径方向成分と接線方向成分のそれぞれに起因する振動値に基づいて、モータの振動現象を算出することが可能となる。これにより、第１の実施形態と同様、従来よりもさらに計算負荷を抑制し、計算の精度を向上させた、モータにおける振動現象を特定するための算出方法を提供することが可能となる。

【0149】

＜その他の実施形態＞
上述した実施形態では、本願発明に係る算出方法の適用事例として、モータの振動現象のシミュレーションを例に挙げて説明した。これに限らず、例えば、モータ制御の際に本願発明に係る算出方法の算出結果を用いてもよい。具体的には、モータにて発生している振動現象（振動値）と、モータの制御のためのパラメータとの対応関係を定義したテーブルを予め規定しておく。ここでのパラメータとしては、モータの回転速度やモータへの負荷に対する最適な電流が規定されていてよい。そして、モータ制御装置の制御部が、本願発明に係る算出方法にて振動値を算出し、その振動値に対応するパラメータを、上記テーブルを参照することで取得してモータの制御を行うような構成であってもよい。

【0150】

本願発明に係る手法をモータの制御に用いる場合には、例えば、回転数制御もしくはトルク制御にて検知される、回転数測定値と電流測定値から、本手法により最適化された上記の２つの測定値に対応する高調波電流の波高値および位相を適用して、より低振動となるように制御することができる。上述したように、例えば、空隙磁束密度基準による手法では、少ティース数のモータの振動の推定には不向きであり、これに対して、本願発明に係る手法は少ティース数のモータの低振動への制御にも適用することができる。

【0151】

また、本願発明において、上述した１以上の実施形態の機能を実現するためのプログラムやアプリケーションを、ネットワーク又は記憶媒体等を用いてシステム又は装置に供給し、そのシステム又は装置のコンピュータにおける１つ以上のプロセッサがプログラムを読出し実行する処理でも実現可能である。

【0152】

また、１以上の機能を実現する回路（例えば、ＡＳＩＣ（Ａｐｐｌｉｃａｔｉｏｎ Ｓｐｅｃｉｆｉｃ Ｉｎｔｅｇｒａｔｅｄ Ｃｉｒｃｕｉｔ）やＦＰＧＡ（Ｆｉｅｌｄ Ｐｒｏｇｒａｍｍａｂｌｅ Ｇａｔｅ Ａｒｒａｙ））によって実現してもよい。

【0153】

このように、本発明は上記の実施形態に限定されるものではなく、実施形態の各構成を相互に組み合わせることや、明細書の記載、並びに周知の技術に基づいて、当業者が変更、応用することも本発明の予定するところであり、保護を求める範囲に含まれる。

【0154】

以上の通り、本明細書には次の事項が開示されている。
（１） モータ（例えば、２００）の振動現象の物理量を算出する算出方法であって、
前記モータが備えるステータコア（例えば、２０１）の突極（例えば、２０４）に作用する時間依存の加振力を、空間および時間の周波数ごとに変換して前記モータの振動現象の物理量を算出する算出工程を有し、
前記算出工程において、前記ステータコアの突極数に分割した領域（例えば、２０４ａ～２０４ｌ）のそれぞれに作用する値を前記加振力として用いる、ことを特徴とする算出方法。
この構成によれば、従来よりもさらに計算負荷を抑制し、計算の精度を向上させて、モータにおける振動現象を特定することが可能となる。

【0155】

（２） 前記モータの振動現象の物理量は、円環次数ごとかつ時間周波数ごとに分解された振動値の分布（例えば、式（１５））、もしくは、前記振動値の分布の円環次数を前記ステータコアの突極数に分割した領域に対応して逆変換した振動値（例えば、式（１６））である、ことを特徴とする（１）に記載の算出方法。
この構成によれば、従来よりもさらに計算負荷を抑制し、計算の精度を向上させて、モータにおける振動現象としての振動値の分布や振動値を特定することが可能となる。

【0156】

（３） 前記モータの振動現象の物理量は、円環次数ごとかつ時間周波数ごとに分解された加振力の分布（例えば、式（１３））と、当該加振力の分布に対応して分解されたコンプライアンスの分布（例えば、式（１４）、図３）とをアダマール積算をして得られる振動変位の分布（例えば、式（１５））、もしくは、当該振動変位の分布の円環次数を前記ステータコアの突極数に分割した領域に対応して逆変換した時間依存の周波数の変位（例えば、式（１６））であることを特徴とする（１）または（２）に記載の算出方法。
この構成によれば、従来よりもさらに計算負荷を抑制し、計算の精度を向上させて、モータにおける振動現象としての振動変位の分布や周波数の変位を特定することが可能となる。

【0157】

（４） 前記モータの振動現象の物理量は、円環次数ごとかつ時間周波数ごとに分解された加振力の分布（例えば、式（１３））と、当該加振力の分布に対応して分解されたアクセラレンスの分布（例えば、式（１４）、図４）とをアダマール積算をして得られる振動加速度の分布（例えば、式（１５））、もしくは、当該振動加速度の分布の円環次数を前記ステータコアの突極数に分割した領域に対応して逆変換した時間依存の周波数の変位（例えば、式（１６））であることを特徴とする（１）または（２）に記載の算出方法。
この構成によれば、従来よりもさらに計算負荷を抑制し、計算の精度を向上させて、モータにおける振動現象としての振動加速度の分布や周波数の変位を特定することが可能となる。

【0158】

（５） 前記算出工程において、前記ステータコアの突極数に分割した領域のそれぞれに作用する値のうち、径方向成分を使用することを特徴とする（１）～（４）のいずれかに記載の算出方法。
この構成によれば、ステータコアの突極数に分割した領域に対応する加振力の径方向の値を使用することで、従来よりもさらに計算負荷を抑制することが可能となる。

【0159】

（６） 前記算出工程において、前記加振力を径方向成分と接線方向成分に分離し、前記径方向成分と前記接線方向成分それぞれに基づいて振動現象の物理量を算出し、前記径方向成分と前記接線方向成分それぞれに基づいて算出された振動現象の物理量を、空間および時間の周波数ごとに合算することで、前記モータの振動現象の物理量を算出することを特徴とする（１）～（４）のいずれかに記載の算出方法。
この構成によれば、ステータコアへの加振力の径方向成分と、接線方向成分の両方に起因する振動量に基づいて、モータの振動現象の物理量を算出することが可能となる。

【0160】

（７） 前記算出工程において、以下の式（１）が用いられる、ことを特徴とする（１）～（６）のいずれかに記載の算出方法。

【0161】

【数50】

【0162】

この構成によれば、従来よりもさらに計算負荷を抑制し、計算の精度を向上させてモータにおける振動現象を特定することが可能となる。

【0163】

（８） （１）～（７）のいずれかに記載の算出方法を実行する処理工程と、
前記処理工程にて算出された計算結果と、前記モータの速度および負荷の少なくとも一方に対応して予め規定された制御パラメータとに基づき、前記モータに対する印加電流を制御する制御工程と、
を有することを特徴とするモータ制御方法。
この構成によれば、従来よりもさらに計算負荷を抑制し、計算の精度を向上させてモータにおける振動現象を特定し、その結果に応じて精度良くモータの制御を行うことが可能となる。

【0164】

（９） モータの振動現象の物理量を算出する算出装置であって、
前記モータが備えるステータコアの突極に作用する時間依存の加振力を取得する取得手段と、
前記取得手段にて取得した加振力を用いて、空間および時間の周波数ごとに変換して前記モータの振動現象の物理量を算出する算出手段と
を有し、
前記取得手段は、前記ステータコアの突極数に分割した領域のそれぞれに作用する値を前記加振力として取得する、ことを特徴とする算出装置。
この構成によれば、従来よりもさらに計算負荷を抑制し、計算の精度を向上させてモータにおける振動現象を特定することが可能となる。

【0165】

（１０） コンピュータに、
モータが備えるステータコアの突極に作用する時間依存の加振力を、空間および時間の周波数ごとに変換して前記モータの振動現象の物理量を算出する算出工程を実行させ、
前記算出工程において、前記ステータコアの突極数に分割した領域のそれぞれに作用する値を前記加振力として用いる、ことを特徴とするプログラム。
この構成によれば、従来よりもさらに計算負荷を抑制し、計算の精度を向上させてモータにおける振動現象を特定することが可能となる。

【符号の説明】

【0166】

１…情報処理装置
１０…ＣＰＵ（Ｃｅｎｔｒａｌ Ｐｒｏｃｅｓｓｉｎｇ Ｕｎｉｔ）
１１…ＲＯＭ（Ｒｅａｄ Ｏｎｌｙ Ｍｅｍｏｒｙ）
１２…ＲＡＭ（Ｒａｎｄｏｍ Ａｃｃｅｓｓ Ｍｅｍｏｒｙ）
１３…ＨＤＤ（Ｈａｒｄ Ｄｉｓｋ Ｄｒｉｖｅ）
１４…入力装置
１５…表示装置
１６…通信装置
２００…モータ
２０１…ステータ
２０２…ロータ
２０３…永久磁石
２０４…ティース

【図1】

【図2】

【図3】

【図4】

【図5】

【図6】

【図7】

【図8】

【図9】

【図10】

【図11】

【図12】

【図13】

【図14】

【図15】

【図16】

【図17】

【図18】

【図19】