(57)【要約】【課題】評価関数値を用いて、制御対象を制御する場合において、高い応答性を確保することができ、商品性を向上させることができる制御装置を提供する。【解決手段】制御装置1のECU2は、複数の線分を互いの間に角度を存する状態で結合した関数値を用いて、走行軌道Tr_skを算出するととも、複数の線分の長さを探索軌道補正ベクトルKθを用いて補正することにより、走行軌道Tr_skを補正し(ステップ4)、走行軌道Tr_skを用いて、走行環境モデルを設定するとともに、自動運転車両3を制御し(ステップ21,22)、走行環境モデルを用いて、評価関数値Jを算出し(ステップ5)、走行軌道Tr_skを変化させたときの評価関数値Jの変化の方向を表す移動平均値Pa_iを算出し(ステップ8)、探索軌道補正ベクトルKθは、移動平均値Pa_iを用いて、走行軌道Tr_skを評価関数値Jがその極値に向かって変化する方向に補正するように算出される(ステップ9)。【選択図】図8
2019028524.tif 000002
【特許請求の範囲】
【請求項１】
複数の線分を互いの間に角度を存する状態で結合した関数値を用いて、制御用パラメータを算出する制御用パラメータ算出手段と、
前記複数の線分の長さ及び前記角度の少なくとも一方を補正値を用いて補正することにより、前記制御用パラメータを補正する補正手段と、
当該補正された制御用パラメータを用いて、制御対象の動作状態、又は当該動作状態及び動作環境を仮想的にモデル化した仮想モデルを設定する仮想モデル設定手段と、
前記補正された前記制御用パラメータを用いて、前記制御対象を制御する制御手段と、
前記仮想モデルを用いて、評価関数値を算出する評価関数値算出手段と、
前記制御用パラメータを変化させたときの前記評価関数値の変化の方向を表す方向値を算出する方向値算出手段と、
を備え、
前記補正手段は、前記方向値を用いて、前記制御用パラメータを前記評価関数値が当該評価関数値の極値に向って変化する方向に補正するように、前記補正値を算出することを特徴とする制御装置。
【請求項２】
前記補正手段は、前記補正値に参照信号値を加えた付加補正値を用いて、前記複数の線分の長さ及び前記角度の少なくとも一方を補正し、
前記方向値算出手段は、前記評価関数値に前記参照信号値を乗算した乗算値を用いて、前記方向値を算出することを特徴とする請求項１に記載の制御装置。
【請求項３】
前記参照信号値は、周期性を有する信号値で構成され、
前記補正手段は、前記複数の線分及び前記角度をそれぞれ補正する複数の前記補正値に対して、互いに異なる周波数の複数の前記参照信号値をそれぞれ加えた複数の前記付加補正値を用いて、前記複数の線分の長さ及び前記角度の少なくとも一方を補正することを特徴とする請求項２に記載の制御装置。
【請求項４】
前記制御対象は、自動運転車両であることを特徴とする請求項１ないし３のいずれかに記載の制御装置。
【請求項５】
前記制御対象は、車両に搭載された原動機であることを特徴とする請求項１ないし３のいずれかに記載の制御装置。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、評価関数値を用いて、制御対象を制御する制御装置に関する。
【背景技術】
【０００２】
従来、評価関数値を用いて、制御対象を制御する制御装置として、特許文献１に記載されたものを本出願人は提案済みであり、この制御装置は、制御対象としての排ガス浄化装置に適用されたものである。この排ガス浄化装置は、尿素ＳＣＲタイプのものであり、エンジンの排気通路に設けられたユリア選択還元触媒と、ユリア選択還元触媒の上流側の排気通路内に尿素水を噴射するユリア噴射弁などを備えている。
【０００３】
この制御装置は、アンモニアセンサで検出された検出アンモニア濃度ＮＨ３ｃｏｎｓに基づいて、ユリア噴射弁によるユリア噴射量Ｇｕｒｅａを算出するものであり、このユリア噴射量Ｇｕｒｅａは、基準噴射量Ｇｕｒｅａ＿ｂｓとＦＢ噴射量Ｇｕｒｅａ＿ｆｂの和として算出される。この制御装置は、同文献の図４に示すように、目標値修正部、同定器、予測器、評価関数値算出器及び極値探索最適化器を備えている。
【０００４】
この目標値修正部では、１次遅れフィルタアルゴリズムなどを用いて、修正目標アンモニア濃度ＤＮＨ３ｃｏｎｓ＿ｔｒｇｔ＿ｍｏｄが算出され、同定器では、逐次型最小２乗法アルゴリズムにより、モデルパラメータａ１，ａ２，ｂ１，ｂ２が算出される。さらに、予測器では、モデルパラメータａ１，ａ２，ｂ１，ｂ２、探索入力ＤＧｕｒｅａ＿ｅｘｓ及び予測アルゴリズムを用いて、予測アンモニア濃度ＰＲＥＤＮＨ３ｅｘｓが算出され、評価関数値算出器では、修正目標アンモニア濃度ＤＮＨ３ｃｏｎｓ＿ｔｒｇｔ＿ｍｏｄと予測アンモニア濃度ＰＲＥＤＮＨ３ｅｘｓとの偏差の二乗和として、評価関数値Ｊが算出される。
【０００５】
また、極値探索最適化器では、負値化した評価関数値−Ｊに、参照信号値Ｓｒｅｆを乗算することにより、両者の積ＣＲを算出し、この積ＣＲの移動平均値ＣＲａｖｅを算出する。そして、移動平均値ＣＲａｖｅが値０になるように、ＦＢ噴射量Ｇｕｒｅａ＿ｆｂに相当する最適噴射量ＤＧｕｒｅａ＿ｏｐｔを算出し、これに参照信号値Ｓｒｅｆを加算することにより、探索入力ＤＧｕｒｅａ＿ｅｘｓを算出するとともに、この探索入力ＤＧｕｒｅａ＿ｅｘｓが予測器に出力される。
【先行技術文献】
【特許文献】
【０００６】
【特許文献１】特許第５０８７１４０号公報
【発明の概要】
【発明が解決しようとする課題】
【０００７】
上記従来の制御装置によれば、逐次型最小２乗法アルゴリズムを用いてモデルパラメータａ１，ａ２，ｂ１，ｂ２を算出し、モデルパラメータａ１，ａ２，ｂ１，ｂ２及び予測アルゴリズムを用いて、予測アンモニア濃度ＰＲＥＤＮＨ３ｅｘｓを算出する必要がある関係上、評価関数値Ｊの算出に時間を要するという問題がある。さらに、極値探索最適化器において、例えば複数の値を算出する場合には、複数の参照信号値を用いる必要があるとともに、参照信号値の変動成分がこれらの複数の値に残留しないようにするためには、移動平均フィルタのタップ数を、全ての参照信号値の最小公倍数、又はその倍数とする必要があり、その結果、タップ数が膨大な値となってしまう。その場合、移動平均フィルタの位相遅れが長大となることで、極値探索最適化器における最適化処理の応答性が大幅に低下してしまうことになる。以上の理由により、高い応答性が要求される制御（例えば、自動運転車両の走行軌道の探索制御）には適用することができなくなり、商品性が低下してしまう。
【０００８】
本発明は、上記課題を解決するためになされたもので、評価関数値を用いて、制御対象を制御する場合において、高い応答性を確保することができ、商品性を向上させることができる制御装置を提供することを目的とする。
【課題を解決するための手段】
【０００９】
上記目的を達成するために、請求項１に係る制御装置１，１Ａは、複数の線分を互いの間に角度を存する状態で結合した関数値（図１８，３２）を用いて、制御用パラメータ（走行軌道Ｔｒ＿ｓｋ、エンジン運転条件関数値ＥＮＬ）を算出する制御用パラメータ算出手段（ＥＣＵ２、走行軌道算出部３０、機関効率算出部６２、エミッション排出量算出部６３、ステップ４）と、複数の線分の長さ及び角度の少なくとも一方を補正値（探索軌道補正ベクトルＫθ、運転条件補正ベクトルＫθ”）を用いて補正することにより、制御用パラメータ（走行軌道Ｔｒ＿ｓｋ、エンジン運転条件関数値ＥＮＬ）を補正する補正手段（ＥＣＵ２、走行軌道算出部３０、機関効率算出部６２、エミッション排出量算出部６３ステップ４）と、補正された制御用パラメータ（走行軌道Ｔｒ＿ｓｋ、エンジン運転条件関数値ＥＮＬ）を用いて、制御対象の動作状態、又は動作状態及び動作環境を仮想的にモデル化した仮想モデル（図１９，３２〜３４）を設定する仮想モデル設定手段（ＥＣＵ２、走行軌道算出部３０、機関効率算出部６２、エミッション排出量算出部６３）と、補正された制御用パラメータ（走行軌道Ｔｒ＿ｓｋ、エンジン運転条件関数値ＥＮＬ）を用いて、制御対象（車両３、原動機５）を制御する制御手段（ＥＣＵ２、ステップ２１〜２２）と、仮想モデルを用いて、評価関数値Ｊ，Ｊ”を算出する評価関数値算出手段（ＥＣＵ２、走行環境モデル推定部４０、評価関数値算出部６４、ステップ５）と、制御用パラメータを変化させたときの評価関数値Ｊ，Ｊ”の変化の方向を表す方向値（移動平均値Ｐａ＿ｉ，Ｐａ＿ｉ”）を算出する方向値算出手段（ＥＣＵ２、極値探索コントローラ５０，７０、ステップ８）と、を備え、補正手段は、方向値（移動平均値Ｐａ＿ｉ，Ｐａ＿ｉ”）を用いて、制御用パラメータ（走行軌道Ｔｒ＿ｓｋ、エンジン運転条件関数値ＥＮＬ）を評価関数値Ｊ，Ｊ”が評価関数値Ｊ，Ｊ”の極値に向って変化する方向に補正するように、補正値（探索軌道補正ベクトルＫθ、運転条件補正ベクトルＫθ”）を算出する（ステップ９）ことを特徴とする。
【００１０】
この制御装置によれば、補正された制御用パラメータを用いて、制御対象の動作状態、又は動作状態及び動作環境を仮想的にモデル化した仮想モデルが設定され、仮想モデルを用いて、評価関数値が算出される。この制御用パラメータは、複数の線分を互いの間に角度を存する状態で結合した関数値を用いて算出され、そのような制御用パラメータの補正は、複数の線分の長さ及び角度の少なくとも一方を補正値を用いて補正することによって実行される。そのため、逐次型最小２乗法アルゴリズム及び予測アルゴリズムを用いて、評価関数値を算出する特許文献１の手法を適用した場合には、多数の値を探索し、それを近似関数に置換して制御用パラメータを算出する必要があるのに対して、これと異なり、探索する値の数を大幅に削減することができ、評価関数値を迅速に算出することができる。
【００１１】
さらに、制御用パラメータを変化させたときの評価関数値の変化の方向を表す方向値が算出され、方向値を用いて、制御用パラメータを評価関数値が評価関数値の極値に向かって変化する方向に補正するように、補正値が算出されるので、この補正値により、評価関数値の勾配の情報を用いることなく、評価関数値がその極値になるように、制御用パラメータを補正することができる。そして、そのように補正された制御用パラメータを用いて、制御対象が制御されるので、非線形な応答特性の制御対象や応答特性が時間的に変化するような制御対象（例えば自動運転車両の走行軌道の制御など）を、評価関数値がその極値になるように適切に制御することができる。以上により、高い応答性が要求される制御（例えば、自動運転車両の走行軌道の探索制御）にも適用することができ、高い商品性を確保することができる。
【００１２】
請求項２に係る発明は、請求項１に記載の制御装置１，１Ａにおいて、補正手段は、補正値（探索軌道補正ベクトルＫθ、運転条件補正ベクトルＫθ”）に参照信号値ｗ＿ｉ，ｗ＿ｉ”を加えた付加補正値（最終探索軌道補正ベクトルＫθ＿ｓｋ、最終運転条件補正ベクトルＫθ＿ｓｋ”）を用いて、複数の線分の長さ及び角度の少なくとも一方を補正し、方向値算出手段は、評価関数値Ｊ，Ｊ”に参照信号値ｗ＿ｉ，ｗ＿ｉ”を乗算した乗算値（中間値Ｐｃ＿ｉ，Ｐｃ＿ｉ”）を用いて、方向値（移動平均値Ｐａ＿ｉ，Ｐａ＿ｉ”）を算出することを特徴とする。
【００１３】
この制御装置によれば、補正値に参照信号値を加えた付加補正値を用いて、複数の線分の長さ及び角度の少なくとも一方を補正するので、制御用パラメータは参照信号値に起因する成分を含んだものとなる。さらに、そのような制御用パラメータを用いて、仮想モデルが設定され、仮想モデルを用いて、評価関数値が算出されるので、評価関数値も参照信号値に起因する成分を含んだものとなる。したがって、価関数値に参照信号値を乗算した乗算値は、両者の相関性を表すものとなるので、そのような乗算値を用いることにより、方向値を精度よく算出することができる。その結果、制御対象を、評価関数値がその極値になるように迅速に制御することができる。
【００１４】
請求項３に係る発明は、請求項２に記載の制御装置１，１Ａにおいて、参照信号値ｗ＿ｉ，ｗ＿ｉ”は、周期性を有する信号値で構成され、補正手段は、複数の線分及び角度をそれぞれ補正する複数の補正値（探索軌道補正ベクトルＫθ、運転条件補正ベクトルＫθ”）に対して、互いに異なる周波数の複数の参照信号値ｗ＿ｉ，ｗ＿ｉ”をそれぞれ加えた複数の付加補正値（最終探索軌道補正ベクトルＫθ＿ｓｋ、最終運転条件補正ベクトルＫθ＿ｓｋ”）を用いて、複数の線分の長さ及び角度の少なくとも一方を補正することを特徴とする。
【００１５】
この制御装置によれば、複数の線分及び角度をそれぞれ補正する複数の補正値に対して、互いに異なる周波数の複数の参照信号値をそれぞれ加えた複数の付加補正値を用いて、複数の線分の長さ及び角度の少なくとも一方が補正されるので、多数の折れ点ポイントを有する関数値を用いて、制御用パラメータを算出することができる。その結果、そのような制御用パラメータを用いて制御対象を制御することによって、制御性をさらに向上させることができる。
【００１６】
請求項４に係る発明は、請求項１ないし３のいずれかに記載の制御装置１において、制御対象は、自動運転車両３であることを特徴とする。
【００１７】
この制御装置によれば、自動運転車両を制御する場合において、請求項１〜３に係る発明の作用効果を得ることができる。
【００１８】
請求項５に係る発明は、請求項１ないし３のいずれかに記載の制御装置１Ａにおいて、制御対象は、車両３に搭載された原動機５であることを特徴とする。
【００１９】
この制御装置によれば、車両に搭載された原動機を制御する場合において、請求項１〜３に係る発明の作用効果を得ることができる。
【図面の簡単な説明】
【００２０】
【図１】本発明の第１実施形態に係る制御装置及びこれを適用した自動運転車両の構成を模式的に示す図である。
【図２】第１実施形態の制御装置の機能的な構成を示すブロック図である。
【図３】他車両のリスクポテンシャルＰｒｉｓｋの算出結果の一例を示す図である。
【図４】歩行者のリスクポテンシャルＰｒｉｓｋの算出結果の一例を示す図である。
【図５】交通参加者が自車両の走行レーン上に存在しない走行環境下でのリスクポテンシャルＰｒｉｓｋの算出結果の一例を示す図である。
【図６】交通参加者が自車両の走行レーン上に存在しない走行環境下でのベネフィットポテンシャルＰｂｎｆの算出結果の一例を示す図である。
【図７】交通参加者が存在する走行環境下でのベネフィットポテンシャルＰｂｎｆ及びリスクポテンシャルＰｒｉｓｋの算出結果の一例を示す図である。
【図８】走行軌道算出部の構成を示すブロック図である。
【図９】タイプ１の走行軌道算出手法で用いる基準走行軌道ｗｂｓを定義した関数値を示す図である。
【図１０】タイプ１の走行軌道算出手法における、ｋｖ＞１，ｋｔｒ＜１が成立しているときの縦方向速度ｖｌの算出結果の一例を示す図である。
【図１１】タイプ１の走行軌道算出手法における、ｋｖ＞１，ｋｔｒ＜１が成立しているときの走行軌道Ｔｒ＿ｓｋの算出結果の一例を示す図である。
【図１２】タイプ１の走行軌道算出手法で用いる基準走行軌道ｗｂｓを定義した関数値の変形例を示す図である。
【図１３】タイプ２の走行軌道算出手法で用いる基準走行軌道ｗｂｓを定義した関数値を示す図である。
【図１４】タイプ２の走行軌道算出手法における、ｋｖ＞１，ｋｔｒ＿１＝０．５，ｋｔｒ＿２＜１が成立しているときの走行軌道Ｔｒ＿ｓｋの算出結果の一例を示す図である。
【図１５】タイプ２の走行軌道算出手法で用いる基準走行軌道ｗｂｓを定義した関数値の変形例を示す図である。
【図１６】タイプ３の走行軌道算出手法で用いる軌道重み関数値Ｗｔｒ＿ｉを示す図である。
【図１７】タイプ３の走行軌道算出手法における、ｋｖ＝１又はｋｖ＞１が成立しているときの走行軌道Ｔｒ＿ｓｋの算出結果の一例を示す図である。
【図１８】タイプ４の走行軌道算出手法を説明するための図である。
【図１９】走行環境モデルの算出結果の一例を示す図である。
【図２０】ベネフィットポテンシャル補正係数Ｋｂｎｆの算出に用いるマップの一例を示す図である。
【図２１】図１９のＡ−Ａ線に沿う時刻でのリスクポテンシャルＰｒｉｓｋ及びベネフィットポテンシャルＰｂｎｆを示す図である。
【図２２】最も理想的な走行軌道Ｔｒ＿ｓｋの位置を説明するための図である。
【図２３】評価関数値Ｊと最終探索軌道補正ベクトルＫθ＿ｓｋの１つの要素ｋθ＿１＿ｓｋとの関係を説明するための図である。
【図２４】移動平均値Ｐａ＿ｉと最終探索軌道補正ベクトルＫθ＿ｓｋの１つの要素ｋθ＿１＿ｓｋとの関係を説明するための図である。
【図２５】走行軌道探索処理を示すフローチャートである。
【図２６】自動運転制御処理を示すフローチャートである。
【図２７】走行軌道探索処理のシミュレーションを実行した走行環境条件を示す図である。
【図２８】第１実施形態の走行軌道探索処理を実行したときの（ａ）自車両などの横位置、（ｂ）評価関数値Ｊ’及び（ｃ）軌道補正係数ｋｔｒのシミュレーション結果を示す図である。
【図２９】比較のために、ベネフィットポテンシャルＰｂｎｆを省略して走行軌道探索処理を実行したときの（ａ）自車両などの横位置、（ｂ）評価関数値Ｊ’及び（ｃ）軌道補正係数ｋｔｒのシミュレーション結果を示す図である。
【図３０】比較のために、最終探索軌道補正ベクトルＫθ＿ｓｋを省略して走行軌道探索処理を実行したときの（ａ）自車両などの横位置、（ｂ）評価関数値Ｊ’及び（ｃ）軌道補正係数ｋｔｒのシミュレーション結果を示す図である。
【図３１】第２実施形態に係る制御装置の機能的な構成を示すブロック図である。
【図３２】エンジン運転条件関数値ＥＮＬの算出結果の一例を示す図である。
【図３３】機関効率Ｉｔａ＿ｅｎｇの算出に用いるマップの一例を示す図である。
【図３４】エミッション排出量Ｍｅｍの算出に用いるマップの一例を示す図である。
【発明を実施するための形態】
【００２１】
以下、図面を参照しながら、本発明の第１実施形態に係る制御装置について説明する。図１に示すように、この制御装置１は、４輪タイプの自動運転車両３に適用されたものであり、ＥＣＵ２を備えている。なお、以下の説明では、この自動運転車両３（制御対象）を「自車両３」という。
【００２２】
このＥＣＵ２には、状況検出装置４、原動機５及びアクチュエータ６が電気的に接続されている。この状況検出装置４は、カメラ、ミリ波レーダー、レーザーレーダ、ソナー、ＧＰＳ及び各種のセンサなどで構成されており、自車両３の位置及び自車両３の進行方向の周辺状況（交通環境や交通参加者など）を表す周辺状況データＤ＿ｉｎｆｏをＥＣＵ２に出力する。
【００２３】
ＥＣＵ２は、後述するように、この状況検出装置４からの周辺状況データＤ＿ｉｎｆｏに基づいて、自車両３の位置及び自車両３の周辺の交通環境や交通参加者などを認識し、自車両３の未来の走行軌道を決定する。
【００２４】
原動機５は、電気モータ及び内燃機関を組み合わせたタイプいわゆるハイブリッドタイプのものであり、後述するように、自車両３の未来の走行軌道が決定されたときに、自車両３がこの走行軌道で走行するように、ＥＣＵ２によって原動機５の出力が制御される。
【００２５】
また、アクチュエータ６は、制動用アクチュエータ及び操舵用アクチュエータなどで構成されており、後述するように、自車両３の未来の走行軌道が決定されたときに、自車両３がこの走行軌道で走行するように、ＥＣＵ２によってアクチュエータ６の動作が制御される。
【００２６】
一方、ＥＣＵ２は、ＣＰＵ、ＲＡＭ、ＲＯＭ、Ｅ２ＰＲＯＭ、Ｉ／Ｏインターフェース及び各種の電気回路（いずれも図示せず）などからなるマイクロコンピュータで構成されており、上述した状況検出装置４からの周辺状況データＤ＿ｉｎｆｏなどに基づいて、後述するように、走行軌道決定処理などを実行する。なお、本実施形態では、ＥＣＵ２が、制御用パラメータ算出手段、補正手段、仮想モデル設定手段、制御手段、評価関数値算出手段及び方向値算出手段に相当する。
【００２７】
次に、図２を参照しながら、本実施形態の制御装置１の機能的な構成について説明する。この制御装置１は、以下に述べる算出アルゴリズムによって、走行軌道Ｔｒ＿ｓｋを算出し、この走行軌道Ｔｒ＿ｓｋで走行するように、自車両３の走行状態を制御するものである。
【００２８】
同図に示すように、制御装置１は、リスクポテンシャル算出部２０、ベネフィットポテンシャル算出部２１及び走行軌道算出部３０を備えており、これらの要素２０，２１，３０は、具体的にはＥＣＵ２によって構成されている。なお、本実施形態では、走行軌道算出部３０が制御用パラメータ算出手段、補正手段及び仮想モデル設定手段に相当し、走行軌道Ｔｒ＿ｓｋが制御用パラメータに相当する。
【００２９】
まず、リスクポテンシャル算出部２０について説明する。このリスクポテンシャル算出部２０では、前述した周辺状況データＤ＿ｉｎｆｏに基づき、交通参加者（歩行者や車両）の現在から未来における存在確率や、走行してはいけない領域を表すマップ（図示せず）を作成し、走行軌道Ｔｒ＿ｓｋに応じて、このマップを検索することにより、リスクポテンシャルＰｒｉｓｋが算出される。このリスクポテンシャルＰｒｉｓｋは、現在から未来において、自車両３の進行方向の周辺における交通参加者が存在する可能性（確率）がある領域、及び自車両３が走行すべきでない走行不可能域が存在する可能性（確率）がある領域などを表す値であり、具体的には、図３〜５に示すように算出される。
【００３０】
まず、図３は、自車両３の進行方向に他車両７が存在する走行環境下での、他車両７のリスクポテンシャルＰｒｉｓｋの算出結果の一例を表している。同図の横軸の値ｘｒは、自車両３の進行方向に直交する横方向において、自車両３と他車両７との相対位置を表している。この相対位置ｘｒは、他車両７の中心位置を値０として、他車両７の左側が負値に設定され、右側が正値に設定される。同図に示すように、他車両７のリスクポテンシャルＰｒｉｓｋは、他車両７の存在確率が低い領域では値０として算出されるとともに、存在確率が高いほど、より大きい正値として算出される。なお、図３〜５では、リスクポテンシャルＰｒｉｓｋは、便宜上、直線を組み合わせた形で表現されている。
【００３１】
また、図４は、自車両３の進行方向に歩行者８が存在する走行環境下での、歩行者８のリスクポテンシャルＰｒｉｓｋの算出結果の一例を表している。同図に示すように、歩行者８のリスクポテンシャルＰｒｉｓｋは、他車両７の場合と同様に、歩行者８の存在確率が低い領域では値０として算出され、存在確率が高いほど、より大きい正値になるように算出される。
【００３２】
さらに、図５は、自車両３の進行方向に交通参加者が存在しない走行環境下でのリスクポテンシャルＰｒｉｓｋの算出結果の一例を表している。この場合のリスクポテンシャルＰｒｉｓｋは、同図に示す歩道９ａなどのような、自車両３が走行すべきでない走行不可領域における走行を回避すべき確率を表す値である。同図において、横軸の値ｘは、自車両３の幅方向すなわち横方向における自車両３の絶対位置を表しており、この絶対位置ｘは、自車両３の中心位置を値０として、自車両３の左側が負値に設定され、右側が正値に設定されている。この点は、後述する図６においても同様である。
【００３３】
同図に示すように、自車両３の進行方向に交通参加者が存在しない走行環境下でのリスクポテンシャルＰｒｉｓｋは、絶対位置ｘに対するマップ値として算出される。より具体的には、リスクポテンシャルＰｒｉｓｋは、車道９ｂと歩道９ａの境界付近において上昇するとともに、歩道９ａ内では車道９ｂから遠ざかるほど、より大きい正値になるように算出される。
【００３４】
次に、ベネフィットポテンシャル算出部２１について説明する。このベネフィットポテンシャル算出部２１では、周辺状況データＤ＿ｉｎｆｏ及び走行軌道Ｔｒ＿ｓｋ（これまで自車両３が走行してきた走行軌道Ｔｒ＿ｓｋ）に応じて、図示しないマップを検索することにより、ベネフィットポテンシャルＰｂｎｆが算出される。このベネフィットポテンシャルＰｂｎｆは、自車両３が走行するときの、進行方向における理想的な走行領域の存在確率を表す値であり、具体的には、例えば、図６〜７に示すように算出される。
【００３５】
この図６は、自車両３の進行方向に交通参加者が存在しない走行環境下でのベネフィットポテンシャルＰｂｎｆの算出結果の一例を表している。同図に示すように、ベネフィットポテンシャルＰｂｎｆは、負値又は値０になるように算出され、最も理想的な走行位置（☆で示すポイント）で最小値となるとともに、最も理想的な走行位置から横方向に遠ざかるほど、その絶対値が減少するように算出される。なお、図６では、ベネフィットポテンシャルＰｂｎｆは、便宜上、直線を組み合わせた形で表現されている。
【００３６】
また、自車両３の実際の走行中は、例えば、図７に示すような、歩行者８ａ及び他車両７ａなどの交通参加者などが存在する走行環境となる。図７に示す走行環境下の場合、ベネフィットポテンシャルＰｂｎｆ及びリスクポテンシャルＰｒｉｓｋは、図中に示すように算出される。
【００３７】
また、ベネフィットポテンシャルＰｂｎｆを、周辺状況データＤ＿ｉｎｆｏ及び走行軌道Ｔｒ＿ｓｋを入力とし、ベネフィットポテンシャルＰｂｎｆを出力とする深層ニューラルネットワークを用い、熟練ドライバの走行データを学習することによって算出してもよい。さらに、逆強化学習の手法を用いて、熟練ドライバの走行データから熟練ドライバの最適経路を判断する関数値を抽出し、この関数値を用いて、ベネフィットポテンシャルＰｂｎｆを算出してもよい。
【００３８】
次に、図８を参照しながら、前述した走行軌道算出部３０について説明する。この走行軌道算出部３０は、以下に述べるように、走行軌道Ｔｒ＿ｓｋを算出（探索）するものであり、この走行軌道Ｔｒ＿ｓｋは、自車両３が現在から未来にかけて走行すべき軌道に相当する。同図に示すように、走行軌道算出部３０は、走行環境モデル推定部４０及び極値探索コントローラ５０を備えている。
【００３９】
この走行環境モデル推定部４０は、後述するように、周辺状況データＤ＿ｉｎｆｏと、極値探索コントローラ５０からの最終探索軌道補正ベクトルＫθ＿ｓｋとを用いて、走行軌道Ｔｒ＿ｓｋを算出し、評価関数値Ｊを算出するとともに、この評価関数値Ｊを極値探索コントローラ５０に出力する。
【００４０】
また、極値探索コントローラ５０は、走行環境モデル推定部４０から入力された評価関数値Ｊを用いて、上述した最終探索軌道補正ベクトルＫθ＿ｓｋを算出し、これを走行環境モデル推定部４０に出力する。なお、本実施形態では、走行環境モデル推定部４０が評価関数値算出手段に相当し、極値探索コントローラ５０が方向値算出手段に相当する。
【００４１】
次に、上述した走行環境モデル推定部４０について説明する。この走行環境モデル推定部４０では、まず、以下に述べるように、周辺状況データＤ＿ｉｎｆｏに応じて、タイプ１〜４の走行軌道算出手法のいずれかを選択し、選択した走行軌道算出手法によって、走行軌道Ｔｒ＿ｓｋが算出される。なお、以下の説明では、自車両３の進行方向の位置を「縦位置」と呼び、幅方向の位置を「横位置」という。
【００４２】
まず、タイプ１の走行軌道算出手法について説明する。このタイプ１の走行軌道算出手法の場合、走行軌道Ｔｒ＿ｓｋの基準となる基準走行軌道ｗｂｓを、図９に示す関数値のように定義する。この基準走行軌道ｗｂｓは、自車両３が基準縦方向速度ｖｌ＿ｂｓで進行方向に走行したときの、自車両３の相対横位置ｘｆと正規化時刻ｔｎとの関係を定義したもの、すなわち自車両３の進路変更速度を定義したものである。
【００４３】
この場合、縦軸の相対横位置ｘｆは、自車両３の中心位置を値０として、左側の相対位置を正値で、右側の相対位置を負値でそれぞれ表したものである。また、横軸の正規化時刻ｔｎは、自車両３が基準縦方向速度ｖｌ＿ｂｓで進行方向に走行した際の未来時刻を表す値であり、正規化時刻ｔｎ＝０は現在の時刻を表している。以上のように、基準走行軌道ｗｂｓは、基準縦方向速度ｖｌ＿ｂｓとともに定義されている関係上、現在からの未来の経過時間である未来時刻をｔｆとし、基準縦方向速度ｖｌ＿ｂｓで走行したときの現在位置に対する相対縦位置をｙｆとした場合、未来時刻ｔｆにおける相対縦位置ｙｆと縦方向速度ｖｌを表すものとなる。したがって、基準走行軌道ｗｂｓは、図９の横軸の正規化時刻ｔｎを、未来時刻ｔｆ又は相対縦位置ｙｆに置き換えて定義することも可能である。
【００４４】
図９に示す基準走行軌道ｗｂｓと、これを補正するための軌道補正係数ｋｔｒ及び縦方向速度補正係数ｋｖを用いて、任意の未来時刻ｔｆ*における走行軌道Ｔｒ＿ｓｋは、以下のように算出される。この軌道補正係数ｋｔｒは、前述した周辺状況データＤ＿ｉｎｆｏに応じて、正負の所定範囲（例えば、ｗｂｓ・ｋｔｒの絶対値が走路幅を超えない範囲）内の値に設定され、基準縦方向速度ｖｌ＿ｂｓが法定速度又はドライバの要求設定速度のどちらか小さい方に設定されたとした場合、縦方向速度補正係数ｋｖは、前述した周辺状況データＤ＿ｉｎｆｏに応じて、１≦ｋｖが成立する範囲内の値に設定される。
【００４５】
まず、下式（１）により、任意の未来時刻ｔｆ*における縦方向速度ｖｌを算出する。
【数１】
2019028524.tif 000003
【００４６】
次いで、下式（２）により、任意の未来時刻ｔｆ*における相対縦位置ｙｆ*を算出する。
【数２】
2019028524.tif 000004
【００４７】
次に、下式（３）により、任意の正規化時刻ｔｎ*を算出する。
【数３】
2019028524.tif 000005
【００４８】
さらに、上式（３）で算出した任意の正規化時刻ｔｎ*に応じて、図９を検索することにより、基準走行軌道ｗｂｓを算出する。
【００４９】
そして、最終的に、下式（４）により、走行軌道Ｔｒ＿ｓｋを算出する。
【数４】
2019028524.tif 000006
【００５０】
以上の手法により算出されるので、走行軌道Ｔｒ＿ｓｋは、軌道補正係数ｋｔｒ、縦方向速度補正係数ｋｖ及び任意の未来時刻ｔｆ*を独立変数とする関数Ｆ（ｋｔｒ，ｋｖ，ｔｆ*）として算出されることになる。ここで、探索軌道補正ベクトルＫθを下式（５）に示すように定義した場合、走行軌道Ｔｒ＿ｓｋは、探索軌道補正ベクトルＫθ及び任意の未来時刻ｔｆ*を独立変数とする関数Ｆ（Ｋθ，ｔｆ*）として算出されることになる。
【数５】
2019028524.tif 000007
【００５１】
以上のタイプ１の算出手法の場合、例えば、ｋｖ＝１，ｋｔｒ＝１が成立しているときには、上式（１），（４）を参照すると明らかなように、ｖｌ＝ｖｌ＿ｂｓ，Ｔｒ＿ｓｋ＝ｗｂｓが成立することになるので、走行軌道Ｔｒ＿ｓｋは、図９に示す基準走行軌道ｗｂｓと同一になるように算出されることになる。
【００５２】
一方、例えば、ｋｖ＞１，ｋｔｒ＜１が成立しているときには、上式（１），（４）を参照すると明らかなように、縦方向速度ｖｌは、基準縦方向速度ｖｌ＿ｂｓよりも大きい値になる（図１０参照）とともに、走行軌道Ｔｒ＿ｓｋの相対横位置は、基準走行軌道ｗｂｓよりも小さい値になる。その結果、走行軌道Ｔｒ＿ｓｋは、図１１に示すように算出されることになる。この図１１は、理解の容易化のために、横軸を相対縦位置ｙｆで表したものである。
【００５３】
同図において、ｙｆ'，ｙｆ*は、ｋｖ＝１，ｋｔｒ＝１が成立しているとき及びｋｖ＞１，ｋｔｒ＜１が成立しているときの、任意の正規化時刻ｔｎ*に対応する相対縦位置を表している。同図に示すように、ｋｖ＞１，ｋｔｒ＜１が成立している場合、ｖｌ＞ｖｌ＿ｂｓが成立することで、相対縦位置ｙｆ*は、ｋｖ＝１，ｋｔｒ＝１が成立しているときの相対縦位置ｙｆ'よりも前方に移動した位置になるとともに、相対横位置もｋｖ＝１，ｋｔｒ＝１が成立しているときよりも小さい値になる。
【００５４】
この場合、図１１は、理解の容易化のために、横軸を相対縦位置ｙｆで表したものであるが、前述したように、横軸を正規化時刻ｔｎや未来時刻ｔｆに置き換えて表すことも可能である。そのため、このタイプ１の走行軌道算出手法の場合には、実際の走行軌道Ｔｒ＿ｓｋは、未来時刻ｔｆと関連付けられた状態で算出される。この点は、後述するタイプ２やタイプ３の走行軌道算出手法などにおいても同様である。
【００５５】
以上のタイプ１の走行軌道算出手法は、例えば、前方の走行車両を追い抜く走行環境下（後述する図２７参照）や、車線変更などを実行する走行環境下で、走行軌道Ｔｒ＿ｓｋを算出するのに最適な手法である。
【００５６】
なお、タイプ１の走行軌道算出手法において、基準走行軌道ｗｂｓとして、図９に示す関数値に代えて、図１２に示すような、車両の操舵特性を考慮した関数値を用いてもよい。
【００５７】
次に、前述したタイプ２の走行軌道算出手法について説明する。このタイプ２の走行軌道算出手法の場合、基準軌道Ｗｂｓを図１３に示す関数値のように定義する。この図１３の関数値と前述した図９の関数値を比較すると明らかなように、図１３の関数値の場合、複数の基準走行軌道ｗｂｓが設定されているとともに、これらの複数の基準走行軌道ｗｂｓは、第１軌道設定値ｋｔｒ＿１に応じたマップ値として設定されている。この第１軌道設定値ｋｔｒ＿１は、走行軌道の左右方向の位置変更の際における移動速度を設定（選択）するためのものであり、前述した周辺状況データＤ＿ｉｎｆｏに応じて、正負の所定範囲（例えば、ｗｂｓ・ｋｔｒ＿１の絶対値が走路幅を超えない範囲）内の値に設定されている。
【００５８】
また、このタイプ２の走行軌道算出手法の場合、前述したタイプ１の走行軌道算出手法における軌道補正係数ｋｔｒに相当する値が、第２軌道補正係数ｋｔｒ＿２として設定される。すなわち、基準縦方向速度ｖｌ＿ｂｓが法定速度又はドライバの要求設定速度のどちらか小さい方に設定されたとした場合、第２軌道補正係数ｋｔｒ＿２は、前述した周辺状況データＤ＿ｉｎｆｏに応じて、−１≦ｋｔｒ＿２≦１が成立する範囲内の値に設定される。そして、図１３に示す基準走行軌道ｗｂｓの関数値、第１軌道設定値ｋｔｒ＿１、第２軌道補正係数ｋｔｒ＿２及び縦方向速度補正係数ｋｖを用いて、任意の未来時刻ｔｆ*における走行軌道Ｔｒ＿ｓｋは、以下のように算出される。
【００５９】
まず、下式（６）により、任意の未来時刻ｔｆ*における縦方向速度ｖｌを算出する。
【数６】
2019028524.tif 000008
【００６０】
次いで、下式（７）により、任意の未来時刻ｔｆ*における相対縦位置ｙｆ*を算出する。
【数７】
2019028524.tif 000009
【００６１】
次に、下式（８）により、任意の正規化時刻ｔｎ*を算出する。
【数８】
2019028524.tif 000010
【００６２】
さらに、上式（８）で算出した任意の正規化時刻ｔｎ*及び第１軌道設定値ｋｔｒ＿１に応じて、図１３を検索することにより、基準走行軌道ｗｂｓを算出する。
【００６３】
そして、最終的に、下式（９）により、走行軌道Ｔｒ＿ｓｋを算出する。
【数９】
2019028524.tif 000011
【００６４】
以上の手法により算出されるので、走行軌道Ｔｒ＿ｓｋは、第１軌道設定値ｋｔｒ＿１、第２軌道補正係数ｋｔｒ＿２、縦方向速度補正係数ｋｖ及び任意の未来時刻ｔｆ*を独立変数とする関数Ｆ（ｋｔｒ＿１，ｋｔｒ＿２，ｋｖ，ｔｆ*）として算出されることになる。ここで、探索軌道補正ベクトルＫθを下式（１０）に示すように定義した場合、走行軌道Ｔｒ＿ｓｋは、探索軌道補正ベクトルＫθ及び任意の未来時刻ｔｆ*を独立変数とする関数Ｆ（Ｋθ，ｔｆ*）として算出されることになる。
【数１０】
2019028524.tif 000012
【００６５】
以上のタイプ２の走行軌道算出手法の場合、例えば、ｋｖ＝１，ｋｔｒ＿１＝０，ｋｔｒ＿２＝１が成立しているときには、上式（６），（９）を参照すると明らかなように、ｖｌ＝ｖｌ＿ｂｓ，Ｔｒ＿ｓｋ＝ｗｂｓが成立することになるので、走行軌道Ｔｒ＿ｓｋは、図１３に示すｋｔｒ＿１＝０のときの基準走行軌道ｗｂｓと同一になるように算出されることになる。
【００６６】
一方、例えば、ｋｖ＞１，ｋｔｒ＿１＝０．５，ｋｔｒ＿２＜１が成立しているときには、上式（６），（９）を参照すると明らかなように、縦方向速度ｖｌは、基準縦方向速度ｖｌ＿ｂｓよりも大きい値になる（前述した図１０参照）とともに、走行軌道Ｔｒ＿ｓｋの相対横位置は、ｋｔｒ＿１＝０．５のときの基準走行軌道ｗｂｓよりも小さい値になる。その結果、走行軌道Ｔｒ＿ｓｋは、図１４に示すように算出されることになる。
【００６７】
すなわち、同図に示すように、ｋｖ＞１，ｋｔｒ＿１＝０．５，ｋｔｒ＿２＜１が成立しているときには、その相対縦位置ｙｆ*は、ｋｖ＝１，ｋｔｒ＿１＝０．５，ｋｔｒ＿２＝１が成立しているときの相対縦位置ｙｆ'よりも前方に移動した位置になるとともに、相対横位置もｋｖ＝１，ｋｔｒ＿１＝０．５，ｋｔｒ＿２＝１が成立しているときよりも小さい値になる。
【００６８】
この場合、図１４は、理解の容易化のために、横軸を相対縦位置ｙｆで表したものであるが、このタイプ２の走行軌道算出手法の場合、実際の走行軌道Ｔｒ＿ｓｋは、未来時刻ｔｆと関連付けられた状態で算出される。
【００６９】
以上のタイプ２の走行軌道算出手法は、タイプ１の走行軌道算出手法と同様に、例えば、前方の走行車両を追い抜く走行環境下や、車線変更などを実行する走行環境下で、走行軌道Ｔｒ＿ｓｋを算出するのに最適な手法である。
【００７０】
なお、タイプ２の走行軌道算出手法において、基準走行軌道ｗｂｓとして、図１３に示す関数値に代えて、図１５に示すような、車両の操舵特性を考慮した関数値を用いてもよい。
【００７１】
次に、前述したタイプ３の走行軌道算出手法について説明する。このタイプ３の走行軌道算出手法の場合、まず、図１６に示すように、４つの軌道重み関数値Ｗｔｒ＿ｉ（ｉ＝１〜４）を正規化時刻ｔｎに対して定義する。
【００７２】
同図に示すように、４つの軌道重み関数値Ｗｔｒ＿ｉは、正規化時刻ｔｎの４つの所定値ｔｎ＿１〜４によって規定される領域に対応して設定されている。第１の軌道重み関数値Ｗｔｒ＿１は、０〜ｔｎ＿２の期間で規定される第１領域に対応し、第２の軌道重み関数値Ｗｔｒ＿２は、ｔｎ＿１〜ｔｎ＿３の期間で規定される第２領域に対応し、第３の軌道重み関数値Ｗｔｒ＿３は、ｔｎ＿２〜ｔｎ＿４と規定される第３領域に対応し、第４の軌道重み関数値Ｗｔｒ＿４は、ｔｎ＿３〜未来と規定される第４領域に対応するように設定されている。
【００７３】
また、４つの軌道重み関数値Ｗｔｒ＿ｉの各々は、上述した対応する領域では値１以下の正値にかつそれ以外の領域では値０に設定されており、隣り合う各２つの軌道重み関数値は、互いにオーバーラップするように設定されているとともに、オーバーラップ部分の２つの軌道重み関数値の和が値１になるように設定されている。
【００７４】
さらに、このタイプ３の走行軌道算出手法では、４つの軌道補正係数ｋｔｒ２＿ｉ（ｉ＝１〜４）が用いられる。これらの４つの軌道補正係数ｋｔｒ２＿ｉは、周辺状況データＤ＿ｉｎｆｏに応じて、正負の所定範囲（例えば、Ｗｔｒ＿ｉ・ｋｔｒ２＿ｉの絶対値が走路幅を超えない範囲）内の値に設定される。
【００７５】
そして、図１６に示す４つの軌道重み関数値Ｗｔｒ＿ｉ、縦方向速度補正係数ｋｖ、４つの軌道補正係数ｋｔｒ２＿ｉを用いて、任意の未来時刻ｔｆ*における走行軌道Ｔｒ＿ｓｋは、以下のように算出される。
【００７６】
まず、下式（１１）により、任意の未来時刻ｔｆ*における縦方向速度ｖｌを算出する。
【数１１】
2019028524.tif 000013
【００７７】
次いで、下式（１２）により、任意の未来時刻ｔｆ*における相対縦位置ｙｆ*を算出する。
【数１２】
2019028524.tif 000014
【００７８】
次に、下式（１３）により、任意の正規化時間ｔｎ*を算出する。
【数１３】
2019028524.tif 000015
【００７９】
さらに、上式（１３）で算出した任意の正規化時間ｔｎ*に応じて、図１６を検索することにより、４つの軌道重み関数値Ｗｔｒ＿ｉを算出する。
【００８０】
そして、最終的に、下式（１４）により、走行軌道Ｔｒ＿ｓｋを算出する。
【数１４】
2019028524.tif 000016
【００８１】
このタイプ３の走行軌道算出手法の場合、走行軌道Ｔｒ＿ｓｋは、以上の手法により算出されるので、４つの軌道補正係数ｋｔｒ２＿ｉ、縦方向速度補正係数ｋｖ及び任意の未来時刻ｔｆ*を独立変数とする関数Ｆ（ｋｔｒ２＿ｉ，ｋｖ，ｔｆ*）として算出されることになる。ここで、探索軌道補正ベクトルＫθを下式（１５）に示すように定義した場合、走行軌道Ｔｒ＿ｓｋは、探索軌道補正ベクトルＫθ及び任意の未来時刻ｔｆ*を独立変数とする関数Ｆ（Ｋθ，ｔｆ*）として算出されることになる。なお、この探索軌道補正ベクトルＫθの値は、前述した極値探索コントローラ５０において、後述する手法により算出される。
【数１５】
2019028524.tif 000017
【００８２】
以上のタイプ３の走行軌道算出手法の場合、走行軌道Ｔｒ＿ｓｋは、ｋｖ＝１が成立しているときには、例えば、図１７の実線で示す値として算出され、ｋｖ＞１が成立しているときには、図１７の破線で示す値として算出される。同図において、ｙｆ＿１〜４は、ｋｖ＝１が成立しているときの正規化時間ｔｎ＿１〜４に対応する相対縦位置ｙｆの値である。
【００８３】
同図において、ｙｆ'，ｙｆ*はそれぞれ、ｋｖ＝１及びｋｖ＞１が成立しているときの、任意の正規化時刻ｔｎ*に対応する相対縦位置を表している。同図に示すように、ｋｖ＞１が成立しているときの相対縦位置ｙｆ*は、ｖｌ＞ｖｌ＿ｂｓが成立することで、ｋｖ＝１が成立しているときの相対縦位置ｙｆ'よりも前方に移動した位置になる。
【００８４】
この場合、図１７は、理解の容易化のために、横軸を相対縦位置ｙｆで表したものであるが、このタイプ３の走行軌道算出手法の場合、実際の走行軌道Ｔｒ＿ｓｋは、未来時刻ｔｆと関連付けられた状態で算出される。
【００８５】
以上のタイプ３の走行軌道算出手法は、繁華街を走行する場合や、交通参加者をジグザグに回避しながら走行する場合（後述する図１９参照）などの、複雑な走行環境下で、走行軌道Ｔｒ＿ｓｋを算出するのに最適な手法である。
【００８６】
なお、タイプ３の走行軌道算出手法においては、軌道重み関数値Ｗｔｒ＿ｉとして４つの値Ｗｔｒ＿１〜４を用いたが、軌道重み関数値Ｗｔｒ＿ｉの数はこれに限らず、２〜３又は５つ以上の関数値を用いてもよい。その場合には、軌道補正係数ｋｔｒ２＿ｉも、軌道重み関数値Ｗｔｒ＿ｉと同じ数の値を用いればよい。
【００８７】
次に、前述したタイプ４の走行軌道算出手法について説明する。このタイプ４の走行軌道算出手法の場合、図１８に示すように、自車両３の相対横位置ｘｆを横軸とし、相対縦位置ｙｆを縦軸とする平面において、ｎ（ｎ≧２）個の線分の各２つを互いの間に角度を付けた状態で組み合わせて、走行軌道Ｔｒ＿ｓｋを決定／算出するものである。
【００８８】
具体的には、まず、下式（１６）により、ローカル軌道距離Ｌｉ（ｉ＝１〜ｎ）を算出する。このローカル軌道距離Ｌｉは、ｎ個の線分の各々の長さに相当する。
【数１６】
2019028524.tif 000018
【００８９】
上式（１６）において、ｋｖ＿ｉ（ｉ＝１〜ｎ）は速度補正係数を、ｖ＿ｂｓは基準道なり方向速度を、Δｔはローカル区間走行時間をそれぞれ表している。このローカル区間走行時間Δｔは、自車両３がローカル軌道距離Ｌｉを走行するのに要する時間である。
【００９０】
次いで、下式（１７）により、任意の未来時刻ｔｆ*における相対縦位置ｙｆ*を算出する。
【数１７】
2019028524.tif 000019
上式（１７）のθｉは、図１８に示すように、各線分の自車両３の現時点での中心線に対する角度である。
【００９１】
次に、下式（１８）により、任意の未来時刻ｔｆ*における相対横位置ｘｆ*を算出する。
【数１８】
2019028524.tif 000020
【００９２】
そして、下式（１９）により、走行軌道Ｔｒ＿ｓｋを算出する。
【数１９】
2019028524.tif 000021
【００９３】
さらに、下式（２０）により、道なり方向速度ｖを算出する。この道なり方向速度ｖは、自車両３の走行方向の速度である。
【数２０】
2019028524.tif 000022
【００９４】
なお、図１８は、３個の線分Ｌ１〜Ｌ３を示した例であるが、線分の数はこれに限らず、２個又は４個以上の線分を用いてもよい。
【００９５】
走行軌道Ｔｒ＿ｓｋは、以上の手法により算出されるので、ｎ個の速度補正係数ｋｖ＿ｉ、ｎ個の角度θｉ及び任意の未来時刻ｔｆ*を独立変数とする関数Ｆ（ｋｖ＿ｉ，θｉ，ｔｆ*）として算出されることになる。ここで、探索軌道補正ベクトルＫθを下式（２１）に示すように定義した場合、走行軌道Ｔｒ＿ｓｋは、探索軌道補正ベクトルＫθ及び任意の未来時刻ｔｆ＊を独立変数とする関数Ｆ（Ｋθ，ｔｆ＊）として算出されることになる。なお、この探索軌道補正ベクトルＫθの値は、前述した極値探索コントローラ５０において、後述する手法により算出される。
【数２１】
2019028524.tif 000023
【００９６】
このタイプ４の走行軌道算出手法の場合、探索軌道補正ベクトルＫθの算出式（１７），（１８）において、角度θｉが含まれているので、走行軌道Ｔｒ＿ｓｋを算出する際、隣り合う各２つの線分の間の角度が探索軌道補正ベクトルＫθによって補正（変更）されることになる。
【００９７】
なお、タイプ４の走行軌道算出手法において、探索軌道補正ベクトルＫθをｎ個の速度補正係数ｋｖ＿ｉのみ又はｎ個の角度θｉのみを要素とするベクトルとして構成してもよい。
【００９８】
以上のように、走行環境モデル推定部４０では、周辺状況データＤ＿ｉｎｆｏに基づき、タイプ１〜４の走行軌道算出手法のいずれかを選択し、選択した走行軌道算出手法を用いて、走行軌道Ｔｒ＿ｓｋが算出される。
【００９９】
なお、以上の説明では、タイプ１〜４の走行軌道算出手法において、探索軌道補正ベクトルＫθを用いて、走行軌道Ｔｒ＿ｓｋを算出するように説明したが、走行軌道Ｔｒ＿ｓｋの実際の演算では、後述する理由により、探索軌道補正ベクトルＫθに代えて、極値探索コントローラ５０で算出された最終探索軌道補正ベクトルＫθ＿ｓｋが用いられる。なお、本実施形態では、探索軌道補正ベクトルＫθが補正値に相当し、最終探索軌道補正ベクトルＫθ＿ｓｋが付加補正値に相当する。
【０１００】
さらに、走行環境モデル推定部４０は、走行軌道Ｔｒ＿ｓｋ、リスクポテンシャルＰｒｉｓｋ及びベネフィットポテンシャルＰｂｎｆに基づいて、例えば、図１９に示す走行環境モデル（仮想モデル）を作成するとともに、この図１９の走行環境モデルを用いて、以下に述べる手法により、評価関数値Ｊを所定の算出周期（探索周期）ΔＴｓｋで算出する。
【０１０１】
同図１９は、交通参加者として、歩行者８ｂ、反対車線で停止中の他車両３ｂ及び反対車線を走行中の他車両３ｃが存在する走行環境下での走行環境モデルを平面的に示したものである。
【０１０２】
この走行環境モデル推定部４０では、まず、図１９中の探索インデックスｚｊ（ｊ＝０〜ｍ：ｍは整数）の時刻又は距離で、走行環境モデルを探索（検索）することにより、探索インデックスｚｊに対応するベネフィットポテンシャルＰｂｎｆ（ｚｊ）及びリスクポテンシャルＰｒｉｓｋ（ｚｊ）を算出する。この場合、探索インデックスｚｊとしては、タイプ１〜３の走行軌道算出手法を実行している場合には、未来時刻ｔｆが用いられ、タイプ４の走行軌道算出手法を実行している場合には、道なり距離ｘが用いられる。
【０１０３】
そして、探索したベネフィットポテンシャルＰｂｎｆ（ｚｊ）及びリスクポテンシャルＰｒｉｓｋ（ｚｊ）を用いて、下式（２２）により、評価関数値Ｊを算出する。
【数２２】
2019028524.tif 000024
【０１０４】
上式（２２）において、ｋは算出時刻（探索時刻）であり、Ｋｂｎｆは、ベネフィットポテンシャル補正係数である。このベネフィットポテンシャル補正係数Ｋｂｎｆは、リスクポテンシャルＰｒｉｓｋに応じて、図２０に示すマップを検索することにより算出される。
【０１０５】
同図のＰｒｉｓｋ１〜２は、Ｐｒｉｓｋ１＜Ｐｒｉｓｋ２が成立するように設定されるリスクポテンシャルＰｒｉｓｋの所定値である。このマップでは、ベネフィットポテンシャル補正係数Ｋｂｎｆは、Ｐｒｉｓｋ≦Ｐｒｉｓｋ１の領域では、値１に設定され、Ｐｒｉｓｋ２≦Ｐｒｉｓｋの領域では、値０に設定されているとともに、Ｐｒｉｓｋ１＜Ｐｒｉｓｋ＜Ｐｒｉｓｋ２の領域では、リスクポテンシャルＰｒｉｓｋが大きいほど、より小さい値に設定されている。以上のようなベネフィットポテンシャル補正係数Ｋｂｎｆの設定理由については後述する。また、以上の評価関数値Ｊの算出手法の原理についても後述する。
【０１０６】
次に、前述した極値探索コントローラ５０について説明する。図８に示すように、極値探索コントローラ５０は、ウォッシュアウトフィルタ５１、ｎ個の参照信号発生器５２＿ｉ（ｉ＝１〜ｎ）、ｎ個の乗算器５３＿ｉ（ｉ＝１〜ｎ）、ｎ個の移動平均フィルタ５４＿ｉ、ｎ個の探索コントローラ５５＿ｉ（ｉ＝１〜ｎ）及びベクトル算出部５６を備えている。
【０１０７】
このウォッシュアウトフィルタ５１では、下式（２３）により、フィルタ値Ｐｗが算出される。
【数２３】
2019028524.tif 000025
【０１０８】
上式（２３）に示すように、フィルタ値Ｐｗは、評価関数値の今回値Ｊ（ｋ）と前回値Ｊ（ｋ−１）の差分として算出される。また、ウォッシュアウトフィルタ５１は、評価関数値Ｊに含まれている、後述する参照信号値ｗ＿ｉに起因する周波数成分を通過させるためのものである。この場合、上式（２３）に代えて、後述するｎ個の参照信号値ｗ＿ｉの周波数成分を通過させるバタワースハイパスフィルタアルゴリズムにより、フィルタ値Ｐｗを算出してもよく、ｎ個の参照信号値ｗ＿ｉの各々の周波数成分のみを通過させるｎ個のバンドパスフィルタアルゴリズムにより、ｎ個のフィルタ値Ｐｗを算出するように構成してもよい。
【０１０９】
また、ｎ個の参照信号発生器５２＿ｉからは、ｎ個の参照信号値ｗ＿ｉ（ｉ＝１〜ｎ）がそれぞれ出力される。これらのｎ個の参照信号値ｗ＿ｉは、互いに異なる周期の周期関数値に設定されており、それらの周期は互いに異なるｎ個の値ｍ＿ｉ（ｉ＝１〜ｎ）と算出周期ΔＴｓｋの積ｍ＿ｉ・ΔＴｓｋに設定されている。さらに、周期関数の波形としては、例えば、正弦波、余弦波、三角波、台形波及び矩形波などが用いられる。
【０１１０】
さらに、ｎ個の乗算器５３＿ｉでは、下式（２４）により、ｎ個の中間値Ｐｃ＿ｉ（ｉ＝１〜ｎ）がそれぞれ算出される。
【数２４】
2019028524.tif 000026
【０１１１】
また、ｎ個の移動平均フィルタ５４＿ｉでは、下式（２５）により、ｎ個の移動平均値Ｐａ＿ｉ（方向値）が算出される。
【数２５】
2019028524.tif 000027
【０１１２】
上式（２５）のｍｌｃｍは、上述した値ｍ＿ｉの最大公倍数である。このように、移動平均値Ｐａ＿ｉのサンプリング個数を値ｍｌｃｍ＋１に設定した理由は、移動平均値Ｐａ＿ｉから参照信号値ｗ＿ｉの周波数成分を除去するためである。
【０１１３】
次いで、ｎ個の探索コントローラ５５＿ｉでは、下式（２６），（２７）に示すスライディングモード制御アルゴリズムにより、ｎ個の補正係数ｋθ＿ｉ（ｉ＝１〜ｎ）が算出される。
【０１１４】
【数２６】
2019028524.tif 000028
【数２７】
2019028524.tif 000029
【０１１５】
上式（２６）のσ＿ｉ（ｉ＝１〜ｎ）は切換関数であり、Ｓ＿ｉ（ｉ＝１〜ｎ）は−１＜Ｓ＿ｉ＜０が成立するように設定される応答指定パラメータである。また、式（２７）のＫｓｋ＿ｉ（ｉ＝１〜ｎ）は、所定のゲインである。上式（２６），（２７）を参照すると明らかなように、ｎ個の補正係数ｋθ＿ｉは、適応則入力のみのスライディングモード制御アルゴリズムによって、ｎ個の移動平均値Ｐａ＿ｉを値０に収束させる機能を有するように算出される。
【０１１６】
この場合、ｎ個の補正係数ｋθ＿ｉは、下式（２８）に示す探索軌道補正ベクトルＫθとして表すことができる。
【数２８】
2019028524.tif 000030
【０１１７】
上式（２８）の場合、その右辺の要素は、前述したタイプ１の走行軌道算出手法を用いる場合には前述した式（５）の右辺の要素に、タイプ２の走行軌道算出手法を用いる場合には、前述した式（１０）の右辺の要素に、タイプ３の走行軌道算出手法を用いる場合には、前述した式（１５）の右辺の要素にそれぞれ相当する。さらに、タイプ４の走行軌道算出手法を用いる場合には、以上の所定値ｎに代えて、その２倍値２ｎを用いて式（２２）〜（２８）の算出が実行されるとともに、その場合の式（２８）の右辺が前述した式（２１）の右辺の要素に相当することになる。
【０１１８】
そして、ベクトル算出部５６では、下式（２９），（３０）により、最終探索軌道補正ベクトルＫθ＿ｓｋが算出される。
【０１１９】
【数２９】
2019028524.tif 000031
【数３０】
2019028524.tif 000032
【０１２０】
上式（２９），（３０）と、前述した式（２７），（２８）を比較すると明らかなように、最終探索軌道補正ベクトルＫθ＿ｓｋの場合、その各要素ｋθ＿ｉ＿ｓｋは、探索軌道補正ベクトルＫθの各要素ｋθ＿ｉに参照信号値ｗ＿ｉをそれぞれ加算した値として算出される。したがって、結果的に、走行軌道Ｔｒ＿ｓｋは、探索軌道補正ベクトルＫθを用いて算出されることになる。例えば、タイプ１，２の走行軌道算出手法の場合には、探索軌道補正ベクトルＫθで基準軌道ｗｂｓを補正することによって算出されることになる。
【０１２１】
以上のように、走行環境モデル推定部４０では、前述した算出手法により、評価関数値Ｊが算出され、極値探索コントローラ５０では、前述した算出手法により、最終探索軌道補正ベクトルＫθ＿ｓｋが算出される。次に、これらの評価関数値Ｊ及び最終探索軌道補正ベクトルＫθ＿ｓｋの算出手法の原理について説明する。
【０１２２】
まず、評価関数値Ｊの算出手法の原理について説明する。例えば、タイプ１〜４の走行軌道算出手法のいずれかで走行軌道Ｔｒ＿ｓｋを算出している場合、図１９のＡ−Ａ線で示す探索インデックスの時刻においては、リスクポテンシャルＰｒｉｓｋ及びベネフィットポテンシャルＰｂｎｆは、例えば、図２１に示すものとなる。
【０１２３】
同図に示すように、リスクポテンシャルＰｒｉｓｋは正値で表され、ベネフィットポテンシャルＰｂｎｆは負値で表されるので、最も理想的な走行軌道Ｔｒ＿ｓｋの位置は、図２２に示すような、リスクポテンシャルＰｒｉｓｋとベネフィットポテンシャルＰｂｎｆの和Ｐｒｉｓｋ＋Ｐｂｎｆが極小値（すなわち最小値）になるポイント（☆で示すポイント）となる。
【０１２４】
したがって、走行軌道Ｔｒ＿ｓｋを最も理想的な走行軌道になるように決定する場合、その評価関数値Ｊを、２つの値の和Ｐｒｉｓｋ＋Ｐｂｎｆを用いて定義することが考えられる。しかしながら、ベネフィットポテンシャルＰｂｎｆの領域とリスクポテンシャルＰｒｉｓｋの領域が平面的に見て重なっている条件下（例えば、図１９の探索インデックスｚｍ−ｘの時刻）では、走行軌道Ｔｒ＿ｓｋを２つの値の和Ｐｒｉｓｋ＋Ｐｂｎｆの最小値の位置になるように決定すると、自車両３が交通参加者に接触する可能性が高まってしまう。
【０１２５】
この事象を回避するために、本実施形態の評価関数値Ｊは、その算出式（２２）に示すように、ベネフィットポテンシャル補正係数ＫｂｎｆをベネフィットポテンシャルＰｂｎｆに乗算した値Ｋｂｎｆ・Ｐｂｎｆを用いるとともに、このベネフィットポテンシャル補正係数Ｋｂｎｆを、前述した図２０に示すように、リスクポテンシャルＰｒｉｓｋが所定値Ｐｒｉｓｋ２よりも大きい領域では値０に設定することによって、評価関数値Ｊの算出において、リスクポテンシャルＰｒｉｓｋを優先させるように構成している。
【０１２６】
これにより、ベネフィットポテンシャルＰｂｎｆの領域とリスクポテンシャルＰｒｉｓｋの領域が平面的に見て重なっている条件下では、走行軌道Ｔｒ＿ｓｋがリスクポテンシャルＰｒｉｓｋの大きい領域と交錯するのを回避でき、自車両３が交通参加者に接触するのを回避できることになる。以上の理由により、評価関数値Ｊは、ベネフィットポテンシャル補正係数Ｋｂｎｆを用いて算出される。
【０１２７】
次に、最終探索軌道補正ベクトルＫθ＿ｓｋの算出原理について説明する。まず、走行軌道Ｔｒ＿ｓｋの最適な軌道を算出するには、リスクポテンシャルＰｒｉｓｋが存在する領域との交錯度合いが最小値となるとともに、ベネフィットポテンシャルＰｂｎｆが存在する領域との交錯度合いが最大値となるように、走行軌道Ｔｒ＿ｓｋを算出すればよいことになる。
【０１２８】
言い換えれば、リスクポテンシャルＰｒｉｓｋとベネフィットポテンシャルＰｂｎｆとの和が最小値になるポイントが、その時刻での走行軌道Ｔｒ＿ｓｋの最適な軌道となるので、評価関数値Ｊが極小値（すなわち最小値）となるように、走行軌道Ｔｒ＿ｓｋを決定すればよいことになる。この場合、走行軌道Ｔｒ＿ｓｋは、前述したように、最終探索軌道補正ベクトルＫθ＿ｓｋを用いて算出される関係上、評価関数値Ｊが極小値となるように、最終探索軌道補正ベクトルＫθ＿ｓｋを算出すればよいことになる。
【０１２９】
したがって、本実施形態の場合、評価関数値Ｊが極小値となるように、最終探索軌道補正ベクトルＫθ＿ｓｋを算出するために、以下の原理を用いている。まず、走行軌道Ｔｒ＿ｓｋは、最終探索軌道補正ベクトルＫθ＿ｓｋを用いて算出される関係上、最終探索軌道補正ベクトルＫθ＿ｓｋに含まれる参照信号値ｗ＿ｉの特性（周期関数）に起因して、所定振幅の振動的な挙動を示すことになる。さらに、そのような走行軌道Ｔｒ＿ｓｋを用いて、リスクポテンシャルＰｒｉｓｋ及びベネフィットポテンシャルＰｂｎｆが算出される関係上、評価関数値Ｊも所定振幅の振動的な挙動を示すことになる。
【０１３０】
ここで、最終探索軌道補正ベクトルＫθ＿ｓｋの１つの要素ｋθ＿１＿ｓｋと評価関数値Ｊの関係が図２３に示す曲線として表されると仮定した場合、参照信号値ｗ＿ｉに起因する評価関数値Ｊの振動的な挙動は、図中の矢印Ｙ１又はＹ２に示すように、ある傾きを持った状態となる。一方、前述した移動平均値Ｐａ＿ｉは、評価関数値Ｊのフィルタ値Ｐｗと参照信号値ｗ＿ｉの積の移動平均値であるので、評価関数値Ｊと参照信号値ｗ＿ｉの相関関数に相当する値となる。
【０１３１】
そのため、相関関数に相当する移動平均値Ｐａ＿ｉが正値であれば、評価関数値Ｊの傾きが正値を示し、移動平均値Ｐａ＿ｉが負値であれば、評価関数値Ｊの傾きは負値を示すことになる。これに加えて、移動平均値Ｐａ＿ｉは、前述した式（２５）で算出されることにより、参照信号値ｗ＿ｉの周波数成分が除去した状態で算出される。以上の理由により、移動平均値Ｐａ＿ｉと要素ｋθ＿１＿ｓｋの関係は、例えば、図２４に示すような単調増加の関数として表すことができる。すなわち、移動平均値Ｐａ＿ｉは、走行軌道Ｔｒ＿ｓｋを変更したときに、評価関数値Ｊが変化する方向を表すことになる。
【０１３２】
したがって、評価関数値Ｊが極小値（最小値）になるように、最終探索軌道補正ベクトルＫθ＿ｓｋを算出するには、図２４に示す関数の傾きが値０になるように、移動平均値Ｐａ＿ｉを算出すればよいことになる。すなわち、移動平均値Ｐａ＿ｉが値０に収束するように、フィードバック制御アルゴリズムを用いて、最終探索軌道補正ベクトルＫθ＿ｓｋを算出すればよいことになる。
【０１３３】
以上の理由により、本実施形態の極値探索コントローラ５０では、フィードバック制御アルゴリズムとしてのスライディングモード制御アルゴリズム（２６），（２７）を含む、式（２３）〜（３０）の算出アルゴリズムを用いて、最終探索軌道補正ベクトルＫθ＿ｓｋが算出される。
【０１３４】
このように評価関数値Ｊの極小値を探索する場合、一般的な極値探索アルゴリズムでは、評価関数値Ｊの傾きのデータが必要となるのに対して、本実施形態の算出アルゴリズムでは、評価関数値Ｊの傾きのデータが不要となる。その結果、本実施形態の算出アルゴリズムの場合、一般的な極値探索アルゴリズムと比べて、演算負荷を低減でき、演算精度を向上させることができるという利点がある。同じ理由により、自動運転装置における走行軌道を決定するシステムのような、時間的に変化するシステムに適用した場合、一般的な極値探索アルゴリズムでは、ロバスト性が低下してしまうのに対して、本実施形態の算出アルゴリズムでは、高いロバスト性を確保することができるという利点がある。
【０１３５】
また、前述したように、移動平均値Ｐａ＿ｉは、走行軌道Ｔｒ＿ｓｋを変更したときに、評価関数値Ｊの傾きすなわち評価関数値Ｊが変化する方向を表す値（方向値）であるので、評価関数値Ｊが最小値になるように（すなわち移動平均値Ｐａ＿ｉが値０に収束するように）、最終探索軌道補正ベクトルＫθ＿ｓｋを算出することは、リスクポテンシャルＰｒｉｓｋが存在する領域との交錯度合いが最小値となるとともに、ベネフィットポテンシャルＰｂｎｆが存在する領域との交錯度合いが最大値となるように、最終探索軌道補正ベクトルＫθ＿ｓｋを算出することに相当する。
【０１３６】
さらに、評価関数値Ｊが、ベネフィットポテンシャル補正係数Ｋｂｎｆを用いて算出されるので、リスクポテンシャルＰｒｉｓｋが存在する領域とベネフィットポテンシャルＰｂｎｆが存在する領域とが交錯している場合には、リスクポテンシャルＰｒｉｓｋが大きいほど、走行軌道Ｔｒ＿ｓｋとリスクポテンシャルＰｒｉｓｋが存在する領域との交錯度合いが減少するように、最終探索軌道補正ベクトルＫθ＿ｓｋが算出されることになる。その結果、自車両３が交通参加者に接触するのを回避できることになる。
【０１３７】
次に、図２５を参照しながら、走行軌道探索処理について説明する。この走行軌道探索処理は、前述した算出手法によって、走行軌道Ｔｒ＿ｓｋ、評価関数値Ｊ及び最終探索軌道補正ベクトルＫθ＿ｓｋなどを算出するものであり、ＥＣＵ２によって、前述した所定の算出周期ΔＴｓｋで実行される。なお、以下の説明において算出される各種の値は、ＥＣＵ２のＥ２ＰＲＯＭ内に記憶されるものとする。
【０１３８】
同図に示すように、まず、ステップ１（図では「Ｓ１」と略す。以下同じ）で、状況検出装置５からの周辺状況データＤ＿ｉｎｆｏを読み込む。
【０１３９】
次いで、ステップ２に進み、Ｅ２ＰＲＯＭ内の走行軌道Ｔｒ＿ｓｋ及び周辺状況データＤ＿ｉｎｆｏに応じて、図示しないマップを検索することにより、リスクポテンシャルＰｒｉｓｋを算出する。この場合、リスクポテンシャルＰｒｉｓｋは、探索インデックスｚ０〜ｚｍまでの値が算出され、それらの値がＥ２ＰＲＯＭ内に記憶されている前回値（探索インデックスｚ０〜ｚｍまでの値）に上書きされる。
【０１４０】
次に、ステップ３で、Ｅ２ＰＲＯＭ内の走行軌道Ｔｒ＿ｓｋ及び周辺状況データＤ＿ｉｎｆｏに応じて、図示しないマップを検索することにより、ベネフィットポテンシャルＰｂｎｆを算出する。この場合、ベネフィットポテンシャルＰｂｎｆは、探索インデックスｚ０〜ｚｍまでの値が算出され、それらの値がＥ２ＰＲＯＭ内に記憶されている前回値（探索インデックスｚ０〜ｚｍまでの値）に上書きされる。
【０１４１】
ステップ３に続くステップ４で、走行軌道Ｔｒ＿ｓｋを算出する。具体的には、周辺状況データＤ＿ｉｎｆｏに基づいて、前述したタイプ１〜４の走行軌道算出手法の算出式、すなわち式（１）〜（５）、式（６）〜（１０）、式（１１）〜（１５）及び式（１６）〜（２１）のいずれかを選択し、選択した算出式における探索軌道補正ベクトルＫθの要素を、最終探索軌道補正ベクトルＫθ＿ｓｋの要素に置き換えて、走行軌道Ｔｒ＿ｓｋを算出する。
【０１４２】
この場合、走行軌道Ｔｒ＿ｓｋは、探索インデックスｚ０〜ｚｍまでの値が算出され、それらの値がＥ２ＰＲＯＭ内に記憶されている前回値（探索インデックスｚ０〜ｚｍまでの値）に上書きされる。すなわち、走行軌道Ｔｒ＿ｓｋは、算出周期ΔＴｓｋで逐次、更新されることになる。
【０１４３】
次いで、ステップ５に進み、前述した式（２２）により、評価関数値Ｊを算出した後、ステップ６で、前述した式（２３）により、フィルタ値Ｐｗを算出する。
【０１４４】
次に、ステップ７で、前述した式（２４）により、ｎ個の中間値Ｐｃ＿ｉ（乗算値）を算出した後、ステップ８に進み、前述した式（２５）により、ｎ個の移動平均値Ｐａ＿ｉを算出する。
【０１４５】
次いで、ステップ９に進み、前述した式（２６）〜（２８）により、探索軌道補正ベクトルＫθを算出する。
【０１４６】
次に、ステップ１０で、前述した式（２９）〜（３０）により、最終探索軌道補正ベクトルＫθ＿ｓｋを算出した後、本処理を終了する。
【０１４７】
以上のように、本実施形態の制御装置１では、所定の算出周期ΔＴｓｋで、最終探索軌道補正ベクトルＫθ＿ｓｋが、逐次更新されるとともに、次の算出タイミングで、そのように更新された最終探索軌道補正ベクトルＫθ＿ｓｋを用いて、走行軌道Ｔｒ＿ｓｋが逐次更新されることになる。
【０１４８】
次に、図２６を参照しながら、自動運転制御処理について説明する。この自動運転制御処理は、自車両３を、算出された走行軌道Ｔｒ＿ｓｋで走行するように制御するものであり、ＥＣＵ２によって、前述した所定の算出周期ΔＴｓｋよりも長い所定の制御周期ΔＴａｄで実行される。
【０１４９】
同図に示すように、まず、ステップ２０で、Ｅ２ＰＲＯＭ内に記憶されている走行軌道Ｔｒ＿ｓｋを読み込む。
【０１５０】
次いで、ステップ２１に進み、自車両３が読み込んだ走行軌道Ｔｒ＿ｓｋで走行するように、原動機５を駆動する。
【０１５１】
次に、ステップ２２で、自車両３が読み込んだ走行軌道Ｔｒ＿ｓｋで走行するようにアクチュエータ６を駆動した後、本処理を終了する。
【０１５２】
次に、以上のように構成された本実施形態の制御装置１による走行軌道決定処理のシミュレーション結果について説明する。以下に述べるシミュレーション結果は、図２７に示すような、先行車両７ｄ及び対向車両７ｅが存在する走行環境条件下において、前述したタイプ１の走行軌道算出手法を用い、先行車両７ｄを追い越すときのものである。
【０１５３】
まず、図２８は、本実施形態の制御装置１による走行軌道決定処理のシミュレーション結果（以下「本願結果」という）を示しており、図２９は、比較のために、ベネフィットポテンシャルＰｂｎｆを省略したときの走行軌道決定処理のシミュレーション結果（以下「第１比較結果」という）を示しており、図３０は、比較のために、最終探索軌道補正ベクトルＫθ＿ｓｋを省略したときの走行軌道決定処理のシミュレーション結果（以下「第２比較結果」という）を示している。
【０１５４】
また、図２８〜図３０におけるＪ’は、自車両３の先端位置のみでの評価関数値であり、下式（３１）により算出される。
【数３１】
2019028524.tif 000033
この式（３１）のｚ’は、自車両３の先端位置での探索インデックスである。
【０１５５】
まず、第２比較結果を参照すると、時刻ｔ２１で、評価関数値Ｊ’が所定値Ｊ１’付近まで急増していることが判る。この所定値Ｊ１’は、自車両３が交通参加者に衝突したことなどを表す値であり、この時刻ｔ２１における評価関数値Ｊ’の急増は、自車両３が先行車両７ｄに衝突したことに起因するものである。
【０１５６】
また、時刻ｔ２２で、車線変更を開始した後、時刻ｔ２３で、評価関数値Ｊ’が再度、所定値Ｊ１’付近まで急増していることが判る。これは、自車両３が対向車両７ｅと衝突したことに起因するものである。
【０１５７】
これに対して、本願結果の場合、時刻ｔ１で、追い越し動作を開始し、車線変更を実行した後、時刻ｔ２で、理想的な走行位置に復帰している。その間、第２比較結果のような他車両と衝突を生じる場合と異なり、評価関数値Ｊ’がほとんど上昇しておらず、走行軌道Ｔｒ＿ｓｋを適切に決定できていることが判る。
【０１５８】
さらに、第１比較結果を参照すると、この第１比較結果の場合、時刻ｔ１１で、追い越し動作を開始し、車線変更を実行した後、時刻ｔ１２で、元の走行位置に復帰している。その間、評価関数値Ｊ’は、本願結果と同様にほとんど上昇しておらず、その点では、走行軌道Ｔｒ＿ｓｋを適切に決定できていることが判る。
【０１５９】
しかしながら、第１比較結果の場合、時刻ｔ１２以降、自車両３が理想的な走行位置から離間しているのに対して、本願結果の場合、前述したように、時刻ｔ２以降、自車両３が理想的な走行位置に復帰していることが判る。すなわち、本願結果の場合、ベネフィットポテンシャルＰｂｎｆを用いることで、自車両３を理想的な走行位置に復帰させることができ、走行軌道Ｔｒ＿ｓｋをより適切に決定できていることが判る。
【０１６０】
以上のように、第１実施形態の制御装置１によれば、タイプ４の走行軌道算出手法の場合、最終探索軌道補正ベクトルＫθ＿ｓｋを用いて、ローカル軌道距離Ｌｉ及び角度θｉを補正することにより、走行軌道Ｔｒ＿ｓｋが算出される。さらに、この走行軌道Ｔｒ＿ｓｋを用いて、図１９に示すような走行環境モデルが設定され、この走行環境モデルを用いて、評価関数値Ｊが算出される。この走行軌道Ｔｒ＿ｓｋは、ｎ個のローカル軌道距離Ｌｉ（ｉ＝１〜ｎ）を互いの間に角度を存する状態で結合した関数値を用いて算出され、走行軌道Ｔｒ＿ｓｋの補正は、最終探索軌道補正ベクトルＫθ＿ｓｋを用いて、ローカル軌道距離Ｌｉ及び角度θｉを補正することによって実行されるので、逐次型最小２乗法アルゴリズム及び予測アルゴリズムを用いて、評価関数値を算出する特許文献１の場合と異なり、評価関数値Ｊを迅速に算出することができる。
【０１６１】
また、前述したように、移動平均値Ｐａ＿ｉが、評価関数値Ｊと参照信号値ｗ＿ｉの相関関数に相当する値となることで、走行軌道Ｔｒ＿ｓｋを変更したときに、評価関数値Ｊが変化する方向を表す値として算出され、この移動平均値Ｐａ＿ｉを用いて、走行軌道Ｔｒ＿ｓｋを評価関数値Ｊがその極値に向かって変化する方向に補正するように、最終探索軌道補正ベクトルＫθ＿ｓｋが算出されるので、この最終探索軌道補正ベクトルＫθ＿ｓｋにより、評価関数値Ｊの勾配の情報を用いることなく、評価関数値Ｊがその極値になるように、走行軌道Ｔｒ＿ｓｋを補正することができる。すなわち、評価関数値Ｊを満たすような最適解を迅速に得ることができる。そして、そのように補正された走行軌道Ｔｒ＿ｓｋで走行するように、自車両３が制御されるので、非線形な応答特性の自車両３を、評価関数値Ｊがその極値になるように迅速かつ適切に制御することができる。その結果、高い商品性を確保することができる。
【０１６２】
なお、第１実施形態は、探索軌道補正ベクトルＫθ及び最終探索軌道補正ベクトルＫθ＿ｓｋの算出において、評価関数値Ｊを用いた例であるが、評価関数値の負値−Ｊを用いてもよい。その場合には、評価関数値の負値−Ｊが極大値（最大値）になるように、探索軌道補正ベクトルＫθ及び最終探索軌道補正ベクトルＫθ＿ｓｋを算出すればよい。すなわち、評価関数値がその極値に向かって変化する方向に走行軌道Ｔｒ＿ｓｋを補正するように、探索軌道補正ベクトルＫθ及び最終探索軌道補正ベクトルＫθ＿ｓｋを算出すればよい。
【０１６３】
また、第１実施形態は、本発明の制御装置１を４輪車両に適用した例であるが、本発明の制御装置は、これに限らず、２輪車両、３輪車両及び５輪以上の車両にも適用可能である。
【０１６４】
次に、図３１を参照しながら、本発明の第２実施形態に係る制御装置１Ａについて説明する。この制御装置１Ａは、以下に述べるように、前述した車両３の原動機５（制御対象）における内燃機関（以下「エンジン」という）の出力を制御するものである。
【０１６５】
この制御装置１Ａの場合、機械的及び電気的なハードウェアの構成は前述した制御装置１とほぼ同様に構成されているので、以下、異なる点を中心に説明する。また、第１実施形態と同じ構成に対しては同じ符号を付すとともに、その説明を省略する。
【０１６６】
この制御装置１Ａの場合、前述した状況検出装置４は、カメラ、ミリ波レーダー、レーザーレーダ、ソナー、ＧＰＳ及び各種のセンサに加えて、カーナビゲーションシステムなどを含んでおり、現時点の走行環境の情報データや、カーナビゲーションシステムに記憶された走行環境条件の履歴情報データを表す走行環境情報データＤ＿ｉｎｆｏ２をＥＣＵ２に出力する。この場合、現時点の走行環境の情報データは、車速、勾配及び渋滞状況などに相当し、走行環境の履歴情報データは、車速、勾配及び渋滞状況などに加えて、通勤利用状況や営業利用状況などを含む。
【０１６７】
図３１に示すように、制御装置１Ａは、エンジン運転モデル推定部６０、極値探索コントローラ７０及びエンジン運転条件算出部８０を備えており、これらの要素６０，７０，８０は、具体的には、ＥＣＵ２によって構成されている。
【０１６８】
このエンジン運転モデル推定部６０は、後述するように、走行環境情報データＤ＿ｉｎｆｏ２と、極値探索コントローラ７０からの最終運転条件補正ベクトルＫθ＿ｓｋ”とを用いて、機関効率Ｉｔａ＿ｅｎｇ及びエミッション排出量Ｍｅｍを算出し、評価関数値Ｊ”を算出するとともに、この評価関数値Ｊ”を極値探索コントローラ７０に出力する。この場合、エミッション排出量Ｍｅｍは、エンジンの排気に含まれるＨＣ，ＣＯ及びＮＯｘの量に相当する。
【０１６９】
また、極値探索コントローラ７０（方向値算出手段）は、前述した極値探索コントローラ５０と同様のアルゴリズムにより、エンジン運転モデル推定部６０から入力された評価関数値Ｊ”を用いて、最終運転条件補正ベクトルＫθ＿ｓｋ”及び運転条件補正ベクトルＫθ”を算出し、これらのベクトルＫθ＿ｓｋ”，Ｋθ”をそれぞれ、エンジン運転モデル推定部６０及びエンジン運転条件算出部８０に出力する。なお、本実施形態では、運転条件補正ベクトルＫθ”が補正値に相当し、最終運転条件補正ベクトルＫθ＿ｓｋ”が付加補正値に相当する。
【０１７０】
さらに、エンジン運転条件算出部８０では、後述する手法により、運転条件補正ベクトルＫθ”を用いて、運転条件関数値ＥＮＬを算出する。
【０１７１】
次に、上述したエンジン運転モデル推定部６０について説明する。このエンジン運転モデル推定部６０は、走行環境情報データＤ＿ｉｎｆｏ２及び最終エンジン運転条件補正ベクトルＫθ＿ｓｋ”を用い、評価関数値Ｊ”を算出するものであり、図３１に示すように、要求エンジン出力算出部６１、機関効率算出部６２、エミッション排出量算出部６３及び評価関数値算出部６４を備えている。
【０１７２】
なお、本実施形態では、機関効率算出部６２及びエミッション排出量算出部６３が制御用パラメータ算出手段、補正手段及び仮想モデル設定手段に相当し、評価関数値算出部６４が評価関数値算出手段に相当する。
【０１７３】
この要求エンジン出力算出部６１では、ｍ（ｍは整数）個の要求エンジン出力Ｗｅｎｇ＿ｒｑ＿ｊ（ｊ＝１〜ｍ）がリングバッファ方式で算出される。具体的には、所定の制御周期で、走行環境情報データＤ＿ｉｎｆｏ２に応じて、図示しないマップを検索することにより、最新の要求エンジン出力Ｗｅｎｇ＿ｒｑを算出し、これをｍ番目の値Ｗｅｎｇ＿ｒｑ＿ｍとして設定（更新）すると同時に、前回の算出タイミングで算出されていたｍ−１個の要求エンジン出力Ｗｅｎｇ＿ｒｑ＿ｍ〜Ｗｅｎｇ＿ｒｑ＿２のデータをそれぞれ、要求エンジン出力Ｗｅｎｇ＿ｒｑ＿ｍ−１〜Ｗｅｎｇ＿ｒｑ＿１として設定する。
【０１７４】
なお、これらのｍ個の要求エンジン出力Ｗｅｎｇ＿ｒｑ＿ｊのサンプリング期間は、制御精度の観点から、極値探索コントローラ７０でのベクトルＫθ＿ｓｋ”，Ｋθ”の算出周期よりも十分長い時間に設定されている。
【０１７５】
以上のように算出されたｍ個の要求エンジン出力Ｗｅｎｇ＿ｒｑ＿ｊ（ｊ＝１〜ｍ）は、要求エンジン出力算出部６１から機関効率算出部６２及びエミッション排出量算出部６３に出力される。
【０１７６】
機関効率算出部６２では、ｍ個の要求エンジン出力Ｗｅｎｇ＿ｒｑ＿ｊ及び最終運転条件補正ベクトルＫθ＿ｓｋ”を用いて、以下に述べる手法により、ｍ個の機関効率Ｉｔａ＿ｅｎｇ＿ｊ（ｊ＝１〜ｍ）が算出される。
【０１７７】
まず、下式（３２）により、ｎ（ｎは２以上の整数）個のローカル運転長さＬｉ”（ｉ＝１〜ｎ）を算出する。
【数３２】
2019028524.tif 000034
【０１７８】
上式（３２）において、Ｌｂｓは基準ローカル運転長さであり、エンジン回転数ＮＥ等に応じて、図示しないマップを検索することにより算出される。また、ｋｖ＿ｉ＿ｓｋ”（ｉ＝１〜ｎ）は、後述する運転長さ補正係数ｋｖ＿ｉ”に参照信号値ｗ＿ｉ”が付加された運転長さ補正係数（以下「付加運転長さ補正係数」という）を表しており、この参照信号値ｗ＿ｉ”は、第１実施形態の参照信号値ｗ＿ｉと同様に設定されている。
【０１７９】
次いで、下式（３３）により、エンジン運転条件関数値ＥＮＬ（制御用パラメータ）を算出する。
【数３３】
2019028524.tif 000035
【０１８０】
上式（３３）のθｉ＿ｓｋ”（ｉ＝１〜ｎ）は、後述するローカル運転条件角度θｉ”に参照信号値ｗ＿ｉ”が付加されたローカル運転条件角度（以下「付加ローカル運転条件角度」という）である。また、最終運転条件補正ベクトルＫθ＿ｓｋ”は、下式（３４）に示すように、付加運転長さ補正係数ｋｖ＿ｉ＿ｓｋ”と付加ローカル運転条件角度θｉ＿ｓｋ”を要素とするベクトルとして定義される。
【数３４】
2019028524.tif 000036
【０１８１】
以上により、エンジン運転条件関数値ＥＮＬは、最終運転条件補正ベクトルＫθ＿ｓｋ”及びエンジン回転数ＮＥを独立変数とする関数Ｆ（Ｋθ＿ｓｋ”，ＮＥ）として算出されることになる。
【０１８２】
以上の式（３２）〜（３３）により、エンジン運転条件関数値ＥＮＬは図３２に示すような値として算出される。このエンジン運転条件関数値ＥＮＬは、エンジン負荷Ｅｌｏａｄ（例えば、トルク）に相当するものであるので、図３２の縦軸はエンジン負荷Ｅｌｏａｄに置き換えることが可能である。
【０１８３】
一方、要求エンジン出力Ｗｅｎｇ＿ｒｑと、エンジン負荷Ｅｌｏａｄ及びエンジン回転数ＮＥとの間においては、下式（３５）が成立する。
【数３５】
2019028524.tif 000037
【０１８４】
したがって、要求エンジン出力算出部６１から入力されたｍ個の要求エンジン出力Ｗｅｎｇ＿ｒｑ＿ｊをこの式（３５）の左辺に代入した場合において、エンジン負荷Ｅｌｏａｄ及びエンジン回転数ＮＥの組み合わせが、式（３５）を満たすと同時に、図３２のエンジン運転条件関数値ＥＮＬ（＝Ｅｌｏａｄ）のデータライン上に位置するものを検索すると、ｍ個のエンジン負荷Ｅｌｏａｄ及びエンジン回転数ＮＥの組み合わせが得られることになる。
【０１８５】
そして、これらのｍ個のエンジン負荷Ｅｌｏａｄ及びエンジン回転数ＮＥの組み合わせに応じて、図３３に示すマップを検索することにより、ｍ個の機関効率Ｉｔａ＿ｅｎｇ＿ｊ（ｊ＝１〜ｍ）が算出される。
【０１８６】
一方、エミッション排出量算出部６３では、ｍ個の要求エンジン出力Ｗｅｎｇ＿ｒｑ＿ｊ及び付加運転条件補正ベクトルＫθ＿ｓｋ”を用いて、機関効率算出部６２と同様の手法により、ｍ個のエミッション排出量Ｍｅｍ＿ｊが（ｊ＝１〜ｍ）が算出される。
【０１８７】
すなわち、前述した式（３２）〜（３３）により、エンジン運転条件関数値ＥＮＬを算出し、式（３５）を満たすように、ｍ個のエンジン負荷Ｅｌｏａｄ及びエンジン回転数ＮＥの組み合わせを算出する。
【０１８８】
そして、これらのｍ個のエンジン負荷Ｅｌｏａｄ及びエンジン回転数ＮＥの組み合わせに応じて、図３４に示すマップを検索することにより、ｍ個のエミッション排出量Ｍｅｍ＿ｊ（ｊ＝１〜ｍ）が算出される。
【０１８９】
また、評価関数値算出部６４では、ｍ個の機関効率Ｉｔａ＿ｅｎｇ＿ｊ及びｍ個のエミッション排出量Ｍｅｍ＿ｊを用いて、下式（３６）により、評価関数値Ｊ”を算出する。この式（３７）のＷｊは、Ｗｊ＞０が成立するように設定される重みパラメータである。
【数３６】
2019028524.tif 000038
【０１９０】
なお、以上の評価関数値Ｊ”の算出処理の実行タイミングは、定期的又は走行停止毎に不定期に実行してもよく、エンジンの運転中、連続的に実行してもよい。また、本実施形態では、図３２〜３４の構成が仮想モデルに相当する。
【０１９１】
さらに、前述した極値探索コントローラ７０では、以上のように算出された評価関数値Ｊ”を用いて、最終運転条件補正ベクトルＫθ＿ｓｋ”及び運転条件補正ベクトルＫθ”が算出される。これらのベクトルＫθ＿ｓｋ”，Ｋθ”は、前述した第１実施形態の極値探索コントローラ５０におけるベクトルＫθ＿ｓｋ，Ｋθと同じアルゴリズムにより算出される。
【０１９２】
すなわち、第１実施形態と同様の手法により、フィルタ値Ｐｗ”を算出し、これにｎ個の参照信号値ｗ＿ｉ”を乗算することにより、ｎ個の中間値Ｐｃ＿ｉ（乗算値）を算出する。さらに、ｎ個の中間値Ｐｃ＿ｉに移動平均フィルタ演算を施すことにより、ｎ個の移動平均値Ｐａ＿ｉ”（方向値）を算出し、これが値０になるように、運転条件補正ベクトルＫθ”を算出する。そして、運転条件補正ベクトルＫθ”の各要素に参照信号値ｗ＿ｉ”を加算することにより、最終運転条件補正ベクトルＫθ＿ｓｋ”が算出される。以上により、ベクトルＫθ＿ｓｋ”，Ｋθ”は、評価関数値Ｊ”が極値（極大値又は極小値）になるように、機関効率Ｉｔａ＿ｅｎｇ及びエミッション排出量Ｍｅｍを補正する値として算出される。
【０１９３】
この場合、運転条件補正ベクトルＫθ”は、下式（３７）に示すように定義される。
【数３７】
2019028524.tif 000039
【０１９４】
この式（３７）において、ｋｖ＿ｉ”は、運転長さ補正係数であり、付加運転長さ補正係数ｋｖ＿ｉ＿ｓｋ”から参照信号値ｗ＿ｉ”を減算した値に相当する。また、θｉ”は、ローカル運転条件角度であり、付加ローカル運転条件角度θｉ＿ｓｋ”から参照信号値ｗ＿ｉ”を減算した値に相当する。
【０１９５】
そして、エンジン運転条件算出部８０では、運転条件補正ベクトルＫθ”を用いて、運転条件関数値ＥＮＬが算出される。具体的には、前述した式（３２）〜（３３）において、付加運転長さ補正係数ｋｖ＿ｉ＿ｓｋ”を運転長さ補正係数ｋｖ＿ｉ”に、付加ローカル運転条件角度θｉ＿ｓｋ”をローカル運転条件角度θｉ”にそれぞれ置き換えた式により、運転条件関数値ＥＮＬを算出する。このように、運転条件関数値ＥＮＬは、制御の安定性を確保するために、参照信号値ｗ＿ｉ”に起因する振動成分を含まない値として算出される。
【０１９６】
エンジン運転条件算出部８０において、運転条件関数値ＥＮＬが以上のように算出されると、この運転条件関数値ＥＮＬに基づいて、原動機５におけるエンジンの吸入空気量、燃料噴射量及び点火時期が制御される。
【０１９７】
以上のように、第２実施形態の制御装置１Ａによれば、最終運転条件補正ベクトルＫθ＿ｓｋ”を用いて、エンジン運転条件関数値ＥＮＬが算出され、このエンジン運転条件関数値ＥＮＬを用いて、図３２〜３４に示すような運転モデルが設定されるとともに、これらの運転モデルを用いて、評価関数値Ｊ”が算出される。このエンジン運転条件関数値ＥＮＬは、ｎ個のローカル運転長さＬｉ”を互いの間に角度を存する状態で結合した関数値を用いて算出され、エンジン運転条件関数値ＥＮＬの補正は、ローカル運転長さＬｉ”及び付加ローカル運転条件角度θｉ＿ｓｋ”を補正することによって実行されるので、逐次型最小２乗法アルゴリズム及び予測アルゴリズムを用いて、評価関数値を算出する特許文献１の場合と異なり、評価関数値Ｊ”を迅速に算出することができる。
【０１９８】
また、前述したように、最終運転条件補正ベクトルＫθ＿ｓｋ”が最終探索軌道補正ベクトルＫθ＿ｓｋと同じ手法により算出されるので、評価関数値Ｊ”の勾配の情報を用いることなく、評価関数値Ｊ”がその極値になるように、エンジン運転条件関数値ＥＮＬを補正することができる。そして、そのように補正されたエンジン運転条件関数値ＥＮＬを用いて、原動機５のエンジンの運転状態が制御されるので、非線形な応答特性のエンジンの運転状態を、評価関数値Ｊ”がその極値になるように適切に制御することができる。その結果、高い商品性を確保することができる。
【０１９９】
なお、第２実施形態は、式（３７）に示す運転条件補正ベクトルＫθ”を用いた例であるが、これに代えて、下式（３８）のように定義された運転条件補正ベクトルＫθ”を用いてもよい。その場合には、式（３８）に示す運転条件補正ベクトルＫθ”の各要素に参照信号値ｗ＿ｉ”をそれぞれ加算した値を要素とする最終運転条件補正ベクトルＫθ＿ｓｋ”を用いればよい。
【数３８】
2019028524.tif 000040
【０２００】
また、運転条件補正ベクトルＫθ”をｎ個の運転長さ補正係数ｋｖ＿ｉ”のみ又はｎ個のローカル運転条件角度θｉ”のみを要素とするベクトルとして構成してもよい。
【０２０１】
さらに、第１及び第２実施形態は、本発明の制御装置を自動運転車両及び車両の原動機に適用した例であるが、本発明の制御装置は、これらに限らず、様々な産業機器にも適用可能である。
【符号の説明】
【０２０２】
１  制御装置
１Ａ  制御装置
２  ＥＣＵ（制御用パラメータ算出手段、補正手段、仮想モデル設定手段、制御手段        、評価関数値算出手段、方向値算出手段）
３  自動運転車両（制御対象）
５  原動機（制御対象）
３０  走行軌道算出部（制御用パラメータ算出手段、補正手段、仮想モデル設定手段）
４０  走行環境モデル推定部（評価関数値算出手段）
５０  極値探索コントローラ（方向値算出手段）
６２  機関効率算出部（制御用パラメータ算出手段、補正手段、仮想モデル設定手段）
６３  エミッション排出量算出部（制御用パラメータ算出手段、補正手段、仮想モデル        設定手段）
６４  評価関数値算出部（評価関数値算出手段）
７０  極値探索コントローラ（方向値算出手段）
Ｔｒ＿ｓｋ  走行軌道（制御用パラメータ）
Ｋθ  探索軌道補正ベクトル（補正値、複数の補正値）
Ｋθ＿ｓｋ  最終探索軌道補正ベクトル（付加補正値、複数の付加補正値）
Ｊ  評価関数値
Ｐａ＿ｉ  移動平均値（方向値）
ｗ＿ｉ  参照信号値
Ｐｃ＿ｉ  中間値（乗算値）
ＥＮＬ  エンジン運転条件関数値（制御用パラメータ）
Ｋθ”  運転条件補正ベクトル（補正値）
Ｋθ＿ｓｋ”  最終運転条件補正ベクトル（付加補正値、複数の付加補正値）
Ｊ”  評価関数値
Ｐａ＿ｉ”  移動平均値（方向値）
ｗ＿ｉ”  参照信号値
Ｐｃ＿ｉ”  中間値（乗算値）
【図１】
2019028524.tif 000041
【図２】
2019028524.tif 000042
【図３】
2019028524.tif 000043
【図４】
2019028524.tif 000044
【図５】
2019028524.tif 000045
【図６】
2019028524.tif 000046
【図７】
2019028524.tif 000047
【図８】
2019028524.tif 000048
【図９】
2019028524.tif 000049
【図１０】
2019028524.tif 000050
【図１１】
2019028524.tif 000051
【図１２】
2019028524.tif 000052
【図１３】
2019028524.tif 000053
【図１４】
2019028524.tif 000054
【図１５】
2019028524.tif 000055
【図１６】
2019028524.tif 000056
【図１７】
2019028524.tif 000057
【図１８】
2019028524.tif 000058
【図１９】
2019028524.tif 000059
【図２０】
2019028524.tif 000060
【図２１】
2019028524.tif 000061
【図２２】
2019028524.tif 000062
【図２３】
2019028524.tif 000063
【図２４】
2019028524.tif 000064
【図２５】
2019028524.tif 000065
【図２６】
2019028524.tif 000066
【図２７】
2019028524.tif 000067
【図２８】
2019028524.tif 000068
【図２９】
2019028524.tif 000069
【図３０】
2019028524.tif 000070
【図３１】
2019028524.tif 000071
【図３２】
2019028524.tif 000072
【図３３】
2019028524.tif 000073
【図３４】
2019028524.tif 000074