        //-----------------------------------------------
        //
        //      cont_altnSource_cht_OF10
        //
        //      OF10用のchtの交互計算sourceCode
        //      originalに対し、thermo定義の
        //          「"thermophysicalProperties" → "physicalProperties"」
        //      に修正。codeOptionsを1行追加。
        //
        //  23/10/26    新規作成
        //  24/04/22    scalar2word:"to_string"から"ostringstream"に修正
        //              桁落ちが発生する為（精度向上）。
        //     05/02    scalar2word:有効桁数が設定できる様に修正
        //-----------------------------------------------
        type    coded;
        name    alternateSource_cht;
        libs    ("libutilityFunctionObjects.so");

        codeData
        #{
            label maxWaitTime_;
            label intervalTemp_;
            label writePrecision_;
            word patchName_;
            word tempField_;

            //  label2word
            //    labelをwordに変換
            word label2word(label n)
            {
                char nn[10];
                sprintf(nn, "%d", n);
                return word(nn);
            }

            //  scalar2word
            //    scalarをwordに変換
            word scalar2word(scalar a)
            {
                //std::string st = std::to_string(a);
                std::ostringstream oss;
                oss.precision(writePrecision_);  //有効桁数
                oss << a;
                std::string st = oss.str();
                return word(st);
            }
            
            //  runCommand
            //    shelの実行
            void runCommand(word comm)
            {
                char ccomm[comm.size() + 1];
                strcpy(ccomm, comm.c_str());
                int err = system(ccomm);
                if (err != 0)
                {
                    FatalErrorInFunction
                        << "error: could not run command '" << comm << "' !!"
                        << abort(FatalError);
                }
            }
        #};

        //  1回のみ実行
        codeRead
        #{
        #};

        codeExecute
        #{
            //  runTimeを定義
            Time& runTime( const_cast<Time&>(mesh().time()) );

            //controlDictの内容読み込み
            maxWaitTime_ = readLabel
                (
                    runTime.controlDict().subDict("couplingToFistr")
                        .subDict("jobWaitingTime")
                        .lookup("maxWaitingTime")
                );
            intervalTemp_ = readLabel
                (
                    runTime.controlDict().subDict("couplingToFistr")
                        .subDict("couplingInterval")
                        .lookup("temperature")
                );
            patchName_ = word
                (
                    runTime.controlDict().subDict("couplingToFistr")
                        .lookup("couplingPatch")
                );
            tempField_ = word
                (
                    runTime.controlDict().subDict("couplingToFistr")
                        .subDict("couplingField")
                        .lookup("tempField")
                );
            writePrecision_ = readLabel
                (
                    runTime.controlDict().lookup("writePrecision")
                );

            //変数の定義
            int pNo = Pstream::myProcNo();
            int np = Pstream::nProcs();
            scalar timeValue = runTime.value();
            scalar dT = runTime.deltaTValue();          //deltaT取得
            scalar beforeTime = timeValue - dT;
            label timeIdx = runTime.timeIndex();
            bool isWriteTime = runTime.writeTime();
            bool isTempCoupling = false;
            label patchId = mesh().boundaryMesh().findPatchID(patchName_);
            word procNo = label2word(pNo);              //processNo
            word nProcs = label2word(np);               //全process数
            word tVal = scalar2word(timeValue);         //現在のtime
            word btVal = scalar2word(beforeTime);       //前のtime
            word timeIndex = label2word(timeIdx);       //timeIndex
            const word dir = "coupling_FrontISTR/data/";
            const word tempFile = dir + patchName_ + "_temp_"
                    + timeIndex + "_" + procNo + "tm";
            const word htFile = dir + patchName_ + "_faceHeatFlux_"
                    + tVal + "_" + procNo + "tm";
            word comm;

            //couplingのtimingを確認
            if ( (timeIdx % intervalTemp_ == 0) or
                 (isWriteTime == true) )
            {
                //couplingのtimingを確認
                isTempCoupling = true;
            }

            //couplingのtimingかどうか確認
            if (isTempCoupling == true)
            {
                //fieldを定義
                volScalarField& T = const_cast<volScalarField&>
                    //( mesh().lookupObject<volScalarField>("T") );
                    ( mesh().lookupObject<volScalarField>(tempField_) );
                scalarField& tempVal = refCast<scalarField>
                    ( T.boundaryFieldRef()[patchId] );
                const scalarField tempIntVal = T.boundaryField()[patchId].patchInternalField();
                const polyPatchList& patches = mesh().boundaryMesh();
                const List<label>& patchCellNos = patches[patchId].faceCells();

                //  -------- 熱流速（heatFlux）を出力 ----------------
                // 1/dx を取得
                const List<scalar>& dxInv = 
                    mesh().deltaCoeffs().boundaryField()[patchId];
                // kappa を取得
                //const basicThermo& thermo = mesh().objectRegistry::lookupObject<basicThermo>("thermophysicalProperties");
                const basicThermo& thermo = mesh().objectRegistry::lookupObject<basicThermo>("physicalProperties");
                const volScalarField& volKappa = thermo.kappa();
                List<scalar> kappa(patchCellNos.size(), 0.0);
                forAll (patchCellNos, i)
                {
                    label celli = patchCellNos[i];
                    kappa[i] = volKappa[celli];
                }
                // 熱流束qを出力、保存
                Info<< "(OF) >>> heat flux value is saved to file..." << endl;
                OFstream writeFaceHeatFlux (htFile);
                forAll (tempVal, pi)
                {
                    scalar q = -kappa[pi] * (tempVal[pi] - tempIntVal[pi]) * dxInv[pi];
                    writeFaceHeatFlux << q << endl;
                } 

                //----------- FrontISTRを起動し、温度を計算 ----------
                comm  = "python3 coupling_FrontISTR/python/";
                comm += "contRun_altnSource_cht.py";
                comm += " " + tVal;       //現在のtime
                comm += " " + btVal;      //前のtime
                comm += " " + procNo;     //procNo
                comm += " " + nProcs;     //全proc数
                comm += " " + timeIndex;  //timeIndex(step数)
                //comm += " &";             //裏で起動
                Info<< "(OF) >>> " << comm << endl;
                runCommand(comm);

                //  ------- 温度を読み込み、patchにセット --------
                Info<< "(OF) >>> waiting data of " << tempFile << "..." << endl;
                label flag = 0;
                for ( int count = 0; count < maxWaitTime_; count++)
                {
                    //  file有無確認
                    IFstream readTemp (tempFile);
                    scalar temp = 0.0;
                    if (readTemp.good())
                    {
                        Info<< "(OF) >>> read and set temperature value to patch '" << patchName_ << "'..." << endl;
                        usleep(20000);      //20msのsleep
                        //温度をpatchにセット
                        forAll (tempVal, pi)
                        {
                            readTemp >> temp;
                            tempVal[pi] = temp;
                        }
                        flag = 1;
                        break;
                    }
                    usleep(10000);          //10msのsleep
                }
                if (flag == 0)
                {
                    FatalErrorInFunction
                        << "error: could not read temperature!!"
                        << abort(FatalError);
                }

                //========== OpenFOAMを再開、次の時間の流体温度を計算 ==================
                //  （codedを終了し、戻る）
            
            }
        #};

        codeInclude
        #{
            #include <unistd.h>
            #include "Time.H"
            #include "fvMesh.H"
            #include "volFields.H"
            #include "surfaceFields.H"
            #include "pointFields.H"
            #include "IFstream.H"
            #include "OFstream.H"
            #include "Pstream.H"
            #include "fluidThermo.H"
        #};

        codeOptions
        #{
            -I$(LIB_SRC)/meshTools/lnInclude \
            -I$(LIB_SRC)/thermophysicalModels/basic/lnInclude \
            -I$(LIB_SRC)/transportModels/compressible/lnInclude \
            -I$(LIB_SRC)/transportModels/lnInclude \
            -I$(LIB_SRC)/physicalProperties/lnInclude
        #};
