!       道路財源の一般財源化　　労働供給内生MODEL 基準ケース（総消費需要＋投資需要＋政府の財・サービス購入；金融資産変換）
!                                        BY KYOJI HASHIMOTO 2001.10.9
!       MERRILL  algolism
        REAL,    PARAMETER :: SEIDO=10**8 ! SHUSOKU SEIDO
        INTEGER, PARAMETER :: ZAI=3 !財の数 
        REAL G(0:3,0:3),E(1:3),NEW(0:3),K(1:3)
        INTEGER COUNT,NLABEL,I,J,KK,LJ,S,L(0:3),WA,J1,JM1,INSATU
        REAL MAX,ST,KG,r,W,TR
10      FORMAT(F12.6,F12.6,F12.6)
        LJ=0
        S=0
        ST=30
        INSATU=0
!      ------------MAIN--------------------------
       DO I=1,ZAI
         k(I)=10
       END DO
       CALL INITIAL(G,K,L,NEW)
        CALL EXD(NEW,E,MAX,NLABEL,ST,INSATU)
        L(0)=0
        LJ=0
55      COUNT=COUNT+1
        KG=0
        WA=0
        SHU=0 
!       HANTEI-------------------------------
        SHU=G(0,0)+G(0,1)+G(0,2)+G(0,3)
        IF (SHU.EQ.1) THEN
          WA=0
          DO  I=0,ZAI
             IF(G(0,I).EQ.0) THEN
                WA=L(I)+WA
                IF (L(I).EQ.0) WA=WA+NLABEL
             END IF
          END DO
        END IF        
        IF(WA.EQ.6) THEN
            CALL RSTART(G,K,ST)
            CALL INITIAL(G,K,L,NEW)
            CALL EXD(NEW,E,MAX,NLABEL,ST,INSATU)
            L(0)=0
            LJ=0
            KG=0
!            WRITE(*,*) ST,SEIDO 
         IF (ST.GT.SEIDO) THEN
!            WRITE(*,*) E(1),E(2),E(3)
            INSATU=1 
            CALL EXD(NEW,E,MAX,MAXJ,ST,INSATU)
            STOP           
         END IF 
        END IF
        DO KK=0,ZAI
          IF (NLABEL.EQ.L(KK)) THEN
           J=KK
           L(LJ)=NLABEL
           J1=J+1
           JM1=J-1  
           IF(J.EQ.0) JM1=3
           IF(J.EQ.3) J1=0
           DO  I=0,ZAI
             NEW(I)=G(I,J1)+G(I,JM1)-G(I,J)
             G(I,J)=NEW(I)
           END DO
             IF (NEW(0).EQ.0)  CALL EXD(NEW,E,MAX,NLABEL,ST,INSATU)
             IF (NEW(0).EQ.1) THEN
             DO  I=1,ZAI
                KG=K(I)-G(I,J)
                IF (KG.GT.0) GOTO 209
             END DO
209            NLABEL=I
             END IF
             LJ=KK
             L(J)=0
             GO TO 55
        END IF 
        END DO
        INSATU=1 
!        WRITE(*,*)  E(1),E(2),E(3)
        CALL EXD(NEW,E,MAX,MAXJ,ST,INSATU)
        STOP           
        END
!      ------------------SUB-----------------------------------------
       SUBROUTINE RSTART(G,K,ST) 
       REAL G(0:3,0:3),K(1:3),ST,T(1:3)
       INTEGER J,I
       DO I=1,3
          T(I)=0
       END DO
       DO  J=0,3
         IF (G(0,J).EQ.0) THEN
            DO I=1,3     
             T(I)=T(I)+G(I,J)
            END DO
         END IF
       END DO
       ST=ST*3
       DO  I=1,3
        K(I)=T(I)
       END DO 
       RETURN   
       END
       SUBROUTINE INITIAL(G,K,L,NEW)
       REAL G(0:3,0:3),NEW(0:3),K(1:3),KG
       INTEGER L(0:3),I,J 
       G(0,0)=0
       DO I=1,3
        G(I,0)=K(I)
        G(0,I)=1
       END DO
       DO  I=1,3
         DO  J=1,3
           G(I,J)=G(I,0)
           IF (J.EQ.I) G(I,J)=G(I,0)-1
         END DO
       END DO
       DO 110 I=1,3
       DO 100 J=1,3
         KG=K(I)-G(I,J)
         IF (KG.GT.0) GOTO 105
100    CONTINUE
105    L(I)=I
110    CONTINUE
       DO  I=1,3
          NEW(I)=G(I,0)
       END DO
       RETURN
       END
!      SUBROUTINE EXCESS DEMAND-------------------------
       SUBROUTINE EXD(NEW,E,MAX,MAXJ,ST,INSATU)
       INTEGER, PARAMETER :: MM=10 !家計の数 
       INTEGER, PARAMETER :: JJ=11 !産業の数
       INTEGER, PARAMETER :: ZAI=10 !消費財の数
       REAL RAMDA(1:MM,1:ZAI) !効用関数パラメータ
       REAL KBAR(1:MM),LBAR(1:MM),X(1:MM,1:ZAI)
       REAL LD(1:JJ),KD(1:JJ),LQ(1:JJ),KQ(1:JJ)!要素需要
       REAL FAI(1:JJ),DELTA(1:JJ) !生産関数パラメータ
       REAL GANMA(1:MM),E(1:3),NEW(0:3),q(1:ZAI),P(1:JJ),SUPPLY(1:JJ)
       REAL Q1(1:MM),Q2(1:ZAI)   !現在財、将来財合成価格
       REAL TAX(1:MM),VAT(1:MM),TLTAX(1:MM),TKTAX(1:MM),TF(1:MM)
       REAL MAX,ST,w,r,DW,DR,SIGMA1,LDA,KDA,LSA,KSA
       REAL TTAX,TITAX,TKLTAX,TVAT,TR
       REAL tL,tk,ROADTAX,ROADREV(1:MM),TROADREV,TOKUTEI
       REAL HAIBUN1,UTILITY(1:MM),AIJ,MYU1,MYU2,MYU3
       REAL ALFA(1:MM)!貯蓄のシェアパラメータ
       REAL PH(1:MM),PHD(1:MM),DI(1:MM),H(1:MM)
       REAL C1(1:MM),C2(1:MM),SAVING(1:MM),TOSI(1:10),TOTALX(1:10)
       REAL LABOR(1:MM),KKK,WT1,PHE,EPSI,BUNSI,BUNBO,BETA(1:MM)
       REAL SITA,KINYU(1:10),hokenryo,TSHAHO,tc
       REAL XG(1:10),GC(1:10),GE,gpara!財サービスの政府購入
       real kansetu(1:10),kanseturev(1:10),TKANSETU
       REAL MYU
       INTEGER I,J,MAXJ,INSATU,OS
       DATA LBAR /717.505039,943.1446423,1147.995328,1209.905554,1395.576487,1538.989666,1698.848342,1844.616505,2075.917105,2557.331999/!労働初期保有量
       DATA KBAR /242.6407764,327.2829077,231.5036538,272.9337497,289.8621759,205.2200446,212.9417829,271.300305,341.2414346,526.7116486/!総資本供給
       DATA FAI/2.011297249,1.599710872,2.204250386,1.3808526,1.380852601,1.774826832,1.462019612,1.774826828,1.77482683,1.774826831,1.892741433/!生産関数パラメータ
       DATA DELTA/0.666325567,0.089117481,0.456999796,0.908622149,0.908622149,0.775903797,0.884281417,0.775903798,0.775903798,0.775903798,0.726092008/
       DATA GANMA/0.101390177,0.106540516,0.11733343,0.114126126,0.093296787,0.092286409,0.06946222,0.084374369,0.115732369,0.105457598/!社会保障給付配分パラメータ
!      消費のシェアパラメータ
       DATA RAMDA/0.265131,0.245747,0.243874,0.243625,0.239742,0.232776,0.222052,0.219549,0.202825,0.194413,&
                  0.118868,0.113677,0.095391,0.088581,0.069771,0.066312,0.056058,0.042264,0.041927,0.035056,&
                  0.076270,0.071588,0.065893,0.066700,0.062414,0.060169,0.057604,0.056002,0.052373,0.049700,&
                  0.032848,0.033171,0.034689,0.034195,0.032742,0.034605,0.033461,0.036157,0.037161,0.037159,&
                  0.043740,0.043967,0.048075,0.050792,0.051436,0.054500,0.053153,0.056995,0.057551,0.067398,&
                  0.038186,0.039476,0.035847,0.037776,0.032822,0.031867,0.027512,0.029221,0.027483,0.026201,&
                  0.115816,0.120791,0.118118,0.105737,0.110960,0.125978,0.119336,0.114538,0.122414,0.116537,&
                  0.035804,0.033093,0.044415,0.047483,0.054620,0.057801,0.057001,0.061231,0.055428,0.050856,&
                  0.084110,0.097592,0.096901,0.098635,0.102037,0.104456,0.101993,0.103665,0.102091,0.112044,&
                  0.189227,0.200898,0.216795,0.226477,0.243458,0.231536,0.271830,0.280380,0.300747,0.310637/
!      貯蓄のシェアパラメータ
         DATA ALFA/0.999515667,0.990272398,0.985052604,0.971250555,0.979635582,0.983654196,0.993553371,0.986845904,0.992181163,0.987128593/!貯蓄シェアパラメータ
!      投資の配分パラメータ
         DATA TOSI/0.001947548,0.43793385,0.024382572,0.127504168,0.003245944,0,0.184028499,0,0,0.220957419/
!      労働供給のパラメータ
         DATA BETA/0.987555087,0.990187407,0.992110312,0.991830766,0.992278397,0.993107438,0.99328047,0.993372321,0.99341619,0.994524956/
!      政府購入の配分パラメータ
         DATA GC/0.0140075,0.1974012,0.0506871,0.0293912,0.0000655,0.5881462,0.0578132,0.0011640,0.0000000,0.0613242/
!      間接税実効税率
         DATA KANSETU/0.047756722,0,0,0,0,0,0.432783876,0,0.010840505,0.048503628/!間接税実効税率

!      ------------INITIAL PARAMETER-----------------
       tc=0.05!消費税
       tl=0 !賃金税
       tk=0.2043   !資本税
       ti=0.2   !利子課税税率
       ty=0.095730493 !所得税率
       F=323.0332331  !課税最低限
       hokenryo=0.1  !保険料
       ROADTAX=0.353962673!道路関係間接税率
       HAIBUN1=0.1303127!社会保障給付対税収比
       GPARA=0.5969038 !政府の財サービス購入予算比率
       TOKUTEI=0.650752088!道路特定財源特別会計充当率
       EPSI=0.4
       SIGMA=0.2
       SITA=0.06734235 !利子収入変換パラメータ
!---------------------------------------------------------
       DO J=1,3
         IF(NEW(J).LE.0) THEN
           MAXJ=J
           RETURN
         END IF 
       END DO
        r=NEW(1)/ST
        W=NEW(2)/ST
        TR=NEW(3)/ST
         R=R/W
         TR=TR/W
         W=1
!      ---------------------------------------------------
!       要素需要
        DO  J=1,JJ !商品
            DW=(1-DELTA(J))*W*(1+tl)
            DR=DELTA(J)*r*(1+tk)
          
            LQ(J)=( (DR/DW)**(1-DELTA(J)) )/FAI(J)
            KQ(J)=( (DW/DR)**DELTA(J)     )/FAI(J)
!       生産者価格
            P(J)=(1+tL)*W*LQ(J)+(1+tk)*r*KQ(J)
 !           WRITE(*,*) J,P(J)
       END DO
!       消費者価格
        DO J=1,ZAI
            q(J)=(1+tc+kansetu(j))*P(J) 
        END DO
!       現在財合成価格、将来財合成価格
        DO I=1,10
           Q1(I)=1 
           Q2(I)=1
           DO  J=1,10
               Q1(I)=Q1(I)*(Q(J)/RAMDA(I,J))**RAMDA(I,J)
               Q2(I)=Q2(I)*((Q(J)/(1+(1-ti)*R))/RAMDA(I,J))**RAMDA(I,J)
           END DO
        END DO
!       合成価格の計算        
        DO I=1,10
              PHD(I)=ALFA(I)**SIGMA*Q1(I)**(1-SIGMA)+(1-ALFA(I))**SIGMA*Q2(I)**(1-SIGMA)
              PH(I)=PHD(I)**(1/(1-SIGMA))
!         WRITE(*,*) I,Q(I)
        END DO
        DO  I=1,MM
          TF(I)=GANMA(I)*TR*HAIBUN1    !移転所得
          KINYU(I)=SITA*KBAR(I)        !家計の利子収入
                                       !政府の財サービスの購入
!          WRITE(*,*) I, TF(I),KINYU(I) 
        END DO
!       労働供給量の決定
        DO I=1,MM
           KKK=((1-BETA(I))/BETA(I))**EPSI
           WT1=((1-ty-hokenryo)*W)**EPSI
           PHE=PH(I)**(1-EPSI)
           BUNSI=KKK*PHE*WT1*LBAR(I)-ty*F-(1-ti)*r*KINYU(I)-TF(I)
           BUNBO=(1-ty-hokenryo)*W+KKK*PHE*WT1
           LABOR(I)=BUNSI/BUNBO
           IF (LABOR(I).LT.0) LABOR(I)=0
        END DO
!STOP
        DO I=1,MM
          TAX(I)=ty*(w*LABOR(I)-F)+ti*r*KINYU(I)  !所得税
          IF (LABOR(I).EQ.0) TAX(I)=ti*r*KINYU(I) 
          IF(TAX(I).LT.0) TAX(I)=0
          DI(I)=(1-hokenryo)*W*LABOR(I)+r*KINYU(I)-TAX(I)+TF(I) !家計可処分所得
          H(I)=DI(I)*PHD(I)**(1/(SIGMA-1))
          MYU=(1-EPSI)/EPSI
          UTILITY(I)=( (1-BETA(I))*H(I)**(-MYU)+BETA(I)*(LBAR(I)-LABOR(I))**(-MYU) )**-(1/MYU) 
       END DO
! --------------------------- 消費者の需要----------------------------------
!      家計の現在消費 C1 と将来消費 C2 の需要関数
        DO I=1,10
         C1(I)=ALFA(I)**SIGMA*DI(I)/(Q1(I)**SIGMA*PHD(I))
         C2(I)=(1-ALFA(I))**SIGMA*DI(I)/(Q2(I)**SIGMA*PHD(I))
!         WRITE(*,*) I,C1(I),C2(I)
        END DO
!       貯蓄 
        TSAVE=0
        DO I=1,10
         SAVING(I)=C2(I)*Q2(I)
         TSAVE=TSAVE+(SAVING(I)+(KBAR(I)-KINYU(I))*r)            !家計総貯蓄
        END DO
!       家計の現在財の需要関数
        DO  I=1,10
           DO J=1,10
                X(I,J)=RAMDA(I,J)*Q1(I)*C1(I)/Q(J)
           END DO
        END DO
!       総消費財需要
        DO  J=1,10
            TOTALX(J)=0
            DO I=1,10
               TOTALX(J)=TOTALX(J)+X(I,J)
            END DO
        END DO
!       税収
        TKLTAX=0!資本税、労働税
        DO I=1,MM
          TLTAX(I)=tl*W*LABOR(I)
          TKTAX(I)=tk*r*KBAR(I)
          TKLTAX=TKLTAX+TLTAX(I)+TKTAX(I)
        END DO
!       消費税,所得税
        TVAT=0
        TITAX=0
        TKANSETU=0
        DO I=1,MM
          VAT(I)=0
          KANSETUREV(I)=0
          DO J=1,ZAI
            VAT(I)=VAT(I)+tc*P(J)*X(I,J)
            KANSETUREV(I)=KANSETUREV(I)+kansetu(j)*p(j)*X(I,J)
          END DO
          TVAT=TVAT+VAT(I)
          TITAX=TITAX+TAX(I) 
          TKANSETU=TKANSETU+KANSETUREV(I)
        END DO
!       道路特定財源税収
        TROADREV=0
        DO I=1,MM
         ROADREV(I)=ROADTAX*P(7)*X(I,7)
         TROADREV=TROADREV+ROADREV(I)
        END DO
        TSHAHO=0
!       社会保険料収入
        DO I=1,10
          TSHAHO=TSHAHO+hokenryo*W*LABOR(I)
        END DO 
        TTAX=TITAX+TVAT+(TKANSETU-TROADREV)+TKLTAX+TROADREV*(1-TOKUTEI)+TSHAHO
        SUPPLY(11)=(TR*(1-HAIBUN1-gpara)+TROADREV*TOKUTEI)/P(11) !公共財需要
        GE=gpara*TR !政府の財・サービス購入
!       総需要（総消費需要＋投資需要＋政府の財・サービス購入）
        DO J=1,10
         SUPPLY(J)=TOTALX(J)+TSAVE*TOSI(J)/P(J)+GE*GC(J)/P(J)
        END DO
        DO J=1,JJ !総労働需要、総資本需要
            LD(J)=LQ(J)*SUPPLY(J)
            KD(J)=KQ(J)*SUPPLY(J)
        END DO
        KDA=0
        LDA=0
        DO J=1,JJ
            LDA=LDA+LD(J)
            KDA=KDA+KD(J)
        END DO  
        KSA=0
        LSA=0
        DO I=1,MM !総労働供給、総資本供給
          KSA=KSA+KBAR(I)
          LSA=LSA+LABOR(I)
        END DO  
        E(1)=KDA-KSA
        E(2)=LDA-LSA
        E(3)=TTAX-TR
!       write(*,*) KDA,KSA
!       write(*,*) LDA,LSA
!       write(*,*) ttax,tr
        MAX=0 
        MAXJ=0
        DO  J=1,3
         IF  (E(J).GT.MAX) THEN
              MAX=E(J)
              MAXJ=J
         END IF
        END DO
!       PRINTOUT
        IF (INSATU.GE.1) THEN
          OPEN(UNIT=1,IOSTAT=OS,FILE='基準均衡.CSV',STATUS='NEW')
          WRITE(*,*) '超過需要'
          WRITE(*,*) E(1),E(2),E(3)
          WRITE(*,*) 'w,r,TR'
          WRITE(*,*) W,R,TR
 !所得階級別税負担
          WRITE(1,*) '階級', ',' ,'労働供給',',' ,'消費',',' ,'効用'
        DO I=1,10
          WRITE(*,*) I, LABOR(I),C1(I),UTILITY(I)
          WRITE(1,*) I, ',',LABOR(I),',' ,C1(I),',',UTILITY(I)
        END DO
          WRITE(1,*)  '階級', ',' ,'所得税',',' ,'消費税',',' ,'個別間接税'
        DO I=1,10
          WRITE(*,*)  I,TAX(I),VAT(I),KANSETUREV(I)
          WRITE(1,*)  I,',' ,TAX(I),',' ,VAT(I),',' ,KANSETUREV(I)
        END DO
        END IF
        RETURN     
        end
