Extended logistic type equations in socio-economic growth modelling

Antonio Edoardo Bargellini
Joint work with D. Ritelli and G. Spaletta

Department of Statistical Sciences, University of Bologna, Italy

#WolframTechConf

Abstract

The study models economic growth via a two-parameter extension of the von Bertalanffy equation with seasonality. 
With constant parameters, an exact hypergeometric solution is found.
In the time-varying case, exact bounds and Mathematica simulations validate the results.

Background and motivation/1

• Neoclassical growth theory (Solow - Swan) [7, 8]:

Y = A  (K^α) (L^(1-α)) = A L (k^α)                          (0 < α < 1 )

• Malthusian population dynamics:

L'(t) =  r  L(t)

• Macroeconomic model of Solow [7] :

k'(t) + (ν + δ)  k(t)=  A  ((k(t))^α)                                  (1)
k(0) = (Subscript[k, 0])

where:

• k = K/L ;  K is the capital;  the labour force L follows a Mathus model  L(t)=(Subscript[L, 0]) (𝕖^(νt));
•  A is the technical progress;
• δ is the depreciation of capital;
• ν is the growth rate of the population.

Background and motivation/2

• Limitation:  in the standard model, the labour force with Malthusian growth is unrealistic.
• Literature background: 
      1)  Donghan [10] -> Logistic Law 
      2)  Mingari & Ritelli [2, 3] -> Exact solution in terms of (Subscript[, 2])(Subscript[F, 1]) (Eulerian Hypergeometric function)
      3) Brida & Maldonado [1] -> Exact solution in terms of (Subscript[, 2])(Subscript[F, 1]) with Von Bertalnnfy population dynamics
• Proposal: population growth via extended von Bertalanffy model and other logistic type models
•  Goal: greater realism and flexibility

Euler integral representation of Gauss Hypergeometric Function:

(Subscript[ , 2])(Subscript[F, 1])(a,b,c; z) =(∑_(n=0)^∞)(((Subscript[(a), n])(Subscript[(b), n]))/(Subscript[(c), n]))((z^n)/(n !)) =((Γ(c))/(Γ(c - b)Γ(b))) (∫_( 0)^(  1))((((t )^(b-1))(((1-t))^(c-b-1)))/(((1 - z t))^a))d t

Objectives

• Integrate extended von Bertalanffy [4] model with Solow growth model

• Von Bertalanffy population dynamics :

L'(t)=  r  ((Subscript[L, ∞])- L(t))

where 
• r controls the speed at which L(t)reaches (Subscript[L, ∞])
• (Subscript[L, ∞]) in the carrying capacity.

• Constant parameter von Bertalanffy model

L'(t)=  (Subscript[r, 1])  (Subscript[L, ∞])-  (Subscript[r, 2])L(t)            (2)

• Time - varying parameter von Bertalanffy model

L'(t)=  r(t)  (Subscript[L, ∞])-  (Subscript[r, 2])L(t)          (3)

• Seek exact solutions or analytical bounds for the solution
• Implementation and testing in Mathematica

Population models/1

Von Bertalanffy model with two constant parameters

The solution to (2) is given by:

L(t)  =  ((Subscript[r, 1])/(Subscript[r, 2]))  (Subscript[L, ∞])-  (e^(-(Subscript[r, 2]) t)) (((Subscript[r, 1])/(Subscript[r, 2])) (Subscript[L, ∞]) - (Subscript[L, 0]))

where L(0) = (Subscript[L, 0]) and

ν(t)  =  ((L'(t))/(L(t))) =   (( (Subscript[r, 2 ])((Subscript[L, ∞])(Subscript[r, 1]) - (Subscript[L, 0]) (Subscript[r, 2])))/((Subscript[L, 0]) (Subscript[r, 2]) + (Subscript[L, ∞])(Subscript[r, 1])((e^((Subscript[r, 2]) t))-1)))                 (4)

Inserting (4) into (1) and setting u =(k^(1-α))  , we have :

u'(t) + (1 - α) (ν(t) + δ)  u(t)=   (1 - α) s A 
u(0) = ((Subscript[k, 0])^(1-α))

⟹      u (t) = (e^( H( t ))) (((Subscript[k, 0])^(1-α)) + (1 - α) s A  (∫_( 0)^(  t))(e^(- H( τ )))  d τ)

where  the saving  propensity  s ∈ [0,1]  is  included, and where :

H( t ) = (α - 1)(∫_( 0)^(  t)) ( ν(τ) + δ )  d τ = (α-1) ((δ-(Subscript[r, 2 ]))t + Log ((((Subscript[L, 0]) (Subscript[r, 2]) + (Subscript[L, ∞])(Subscript[r, 1])((e^((Subscript[r, 2]) t))-1))/((Subscript[L, 0]) (Subscript[r, 2])))))

Population models/2

Von Bertalanffy model with two constant parameters

To obtain solution u(t) explicitly, we must compute ( p = (e^((Subscript[r, 2])τ)) ) :

Z( t ) =(∫_( 0)^(  t))(e^(- H( τ )))  d τ =B(∫_( 1)^(  (e^( (Subscript[r, 2])t))))  (((M+ N p ))^(1-α))(p^C) dp=B ( I(t) - I (0))

where

I( t ) =(∫_( 0)^(  (e^( (Subscript[r, 2])t)))) (((M+ N p ))^(1-α))(p^C) dp  =  (((e^((Subscript[r, 2])t(C+1) )) (((M+ N (e^((Subscript[r, 2])t)) ))^(2-α)))/(M (C+1)))    (Subscript[ , 2])(Subscript[F, 1])(1,ℰ+1,ℰ+α ; ((-N (e^((Subscript[r, 2])t)))/M))

and

B = (Subscript[L, 0]^( α-1)) (Subscript[r, 2]^(α-2)),      M = (Subscript[L, 0]) (Subscript[r, 2])-(Subscript[L, ∞ ])(Subscript[r, 1]),      N = (Subscript[L, ∞ ])(Subscript[r, 1]),     ℰ = ((δ (1 - α) )/(Subscript[r, 2])),     C = ℰ + α - 2

⟹      u (t) = (e^(- H( t ))) (((Subscript[k, 0])^(1-α)) + (1 - α)  s  A  B ( I(t) - I (0)) )  
⟹      k(t) = (((u (t)))^(1/(1-α)))

Population models/3

Von Bertalannfy model with time varying parameter

The solution to (3),  with L(0) = (Subscript[L, 0])  and   r(t) = (Subscript[r, 1]) + b Cos(ω t),   is given by:

L(t)  =  ((G(t))/(Subscript[r, 2])) ((e^(-(Subscript[r, 2]) t))/((((Subscript[r, 2])^2) + (ω^2))))                     (5)

G(t) =  b (Subscript[L, ∞])(Subscript[r, 2  ])g(t) (e^((Subscript[r, 2]) t)) -b (Subscript[L, ∞])((Subscript[r, 2])^2) + ((Subscript[L, 0])(Subscript[r, 2])- (Subscript[L, ∞])(Subscript[r, 1])+ (Subscript[L, ∞])(Subscript[r, 1])(e^((Subscript[r, 2]) t))) (((Subscript[r, 2])^2) + (ω^2))

g(t) = (Subscript[r, 2  ]) Cos(ω t)+ω Sin(ω t)

The choice  of  r( t ) follows Essington et al. [11], Katsanevakis & Maravelias [12] , Somers [13].

The growth rate related to (5) is :

ν(t)  =  (( L'(t))/(L(t))) = ((Subscript[r, 2])/(G(t)))(b (Subscript[L, ∞]) ω  h(t)  (e^((Subscript[r, 2]) t)) + b (Subscript[L, ∞])((Subscript[r, 2])^2) - ((Subscript[L, 0])(Subscript[r, 2])- (Subscript[L, ∞])(Subscript[r, 2]))(((Subscript[r, 2])^2) + (ω^2)) )                   (6)

h(t) =  ω Cos(ω t) - (Subscript[r, 2])Sin(ω t)

Population models/4

Von Bertalannfy model with time varying parameter

As before, to linearise (1) with  ν(t) given by  (6),  we use  u =(k^(1-α))  to obtain the Cauchy Problem:

{ u'(t) + (1 - α) (ν(t) + δ)  u(t)=   (1 - α) s A  ,  u(0) = ((Subscript[k, 0])^(1-α))  }

solved by    u (t) = (e^( H( t ))) (((Subscript[k, 0])^(1-α)) + (1 - α) s A  (∫_( 0)^(  t))(e^(- H( τ )))  d τ)  .

Again, we are faced with computing:

Z(t) =(∫_( 0)^(  1))(e^(- H( τ )))  d τ

(e^(- H(τ))) =((((((Subscript[L, 0 ])(Subscript[r, 2]) (((Subscript[r, 2])^2) + (ω^2))(e^(( (Subscript[r, 2])- δ)  τ )) )/(G(τ))) ))^(α-1))

g(t) -> bounded, i.e., admits Max and Min in [0, t] with t > 0 :

(Subscript[g, min])= -(Sqrt[P]),   (Subscript[g, max])= (Sqrt[P]),   P = ((Subscript[r, 2])^2) + (ω^2)

Z( t ) cannot be given explicitly, but we can find explicit  inferior/superior solutions (Subscript[Z, min])(t) and (Subscript[Z, max])(t),   in terms of (Subscript[, 2])(Subscript[F, 1]) ,  bounding Z( t )
•  Change variable p = (e^((Subscript[r, 2])t))
•  Compute (Subscript[Z, min])(t) and (Subscript[Z, max])(t) in terms of (Subscript[I, min])(t) and (Subscript[I, max])(t), similarly to the 2–parameters case
•  Use of Integral representation theorem of  (Subscript[, 2])(Subscript[F, 1])

Implementation and simulations

•  Developed in Mathematica
•  Numerical simulations and parameter sensitivity analysis
•  Comparative plots
•  Validation: improved realism over traditional models

Implementation and simulations 
Constant parameter model

Manipulate[Module[{condizioniOK$},condizioniOK$=r1>r2;
If[condizioniOK$,Plot[{(r1/r2)*Linf-Exp[(-r2)*t]*((r1/r2)*Linf-L0)},{t,initialTime,finalTime},PlotLabels->"L(t)",PlotRange->{0,All},AxesLabel->TraditionalForm/@{t,"L(t)"},LabelStyle->{FontFamily->"Times",FontSize->14},AxesStyle->Arrowheads[{0.,0.03}]],Style[Row[{"⚠️ Conditions not respected: \!\(\*SubscriptBox[\(r\), \(1\)]\)  >  \!\(\*SubscriptBox[\(r\), \(2\)]\) > 0 "," \!\(\*FractionBox[SubscriptBox[\(r\), \(1\)], SubscriptBox[\(r\), \(2\)]]\) = ",r1/r2}],Red,14,Bold]]],{{r1,0.5501},0.0001,1,0.01},{{r2,0.0601},0.0001,1,0.01},{{Linf,5.7601},0.0001,10,0.01},{{L0,3.08},0.1,10,0.01},{{initialTime,0},0},{{finalTime,100.},1,100},Initialization:>{initialTime=0,finalTime=90}]

Implementation and simulations 
Constant parameter model

Manipulate[Module[{condizioniOK,expr},(*controllo vincoli*)condizioniOK=r1 > r2 >0 && L0 r2 !=LInf r1 && (-1+α+(((1-α) δ)/r2)) != 0;
If[condizioniOK,(*definizione espressione*)expr=((((E^(r2 t (-1+α)+t (1-α) δ)) (((L0 r2))^(-1+α)) (((-LInf r1+(E^(r2 t)) LInf r1+L0 r2))^(1-α)) ((k0^(1-α))+A (L0^(-1+α)) (r2^(-2+α)) s (1-α) (-(((((L0 r2))^(2-α)) Hypergeometric2F1[1,1+(((1-α) δ)/r2),α+(((1-α) δ)/r2),-((LInf r1)/(-LInf r1+L0 r2))])/((-LInf r1+L0 r2) (-1+α+(((1-α) δ)/r2))))+((E^(r2 t (-1+α+(((1-α) δ)/r2)))) (((-LInf r1+(E^(r2 t)) LInf r1+L0 r2))^(2-α)) Hypergeometric2F1[1,1+(((1-α) δ)/r2),α+(((1-α) δ)/r2),-(((E^(r2 t)) LInf r1)/(-LInf r1+L0 r2))])/((-LInf r1+L0 r2) (-1+α+(((1-α) δ)/r2)))))))^(1/(1-α)));
(*grafico*)Plot[Evaluate@{expr},{t,initialTime,finalTime},PlotLabels->{"k(t)"},PlotRange->{0,All}, PlotStyle->{Black}, AxesLabel->TraditionalForm/@{t,"k(t)"},LabelStyle->{FontFamily->"Times",FontSize->14},AxesStyle->Arrowheads[{-0.,0.030}]],(*altrimenti:solo messaggio*)Style[Row[{"⚠️ Conditions not respected: (Subscript[r, 1]) > (Subscript[r, 2]) > 0, 0 < δ < 1, 0 < s < 1, 0 < α < 1 and L0 (Subscript[r, 2]) > (Subscript[L, ∞]) (Subscript[r, 1])"," ((Subscript[r, 1])/(Subscript[r, 2])) = ",r1/r2," , δ = ",δ," , s = ",s," , α = ",α, " , L0 (Subscript[r, 2]) - (Subscript[L, ∞]) (Subscript[r, 1]) = ", L0 r2 -LInf r1}],Red,14, Bold]]],{{r1,0.1},0.1,1,0.01},{{r2,0.1},0.1,1,0.01},{{α,0.5},0.1,0.99,0.01},{{L0,0.1},0.1,10,0.01},{{δ,0.1},0.1,0.99,0.01},{{k0,0.5},0.1,10,0.01},{{LInf,1},0.1,10,0.1},{{A,0.5},0.1,10,0.01},{{s,0.5},0.1,0.99,0.01},{{initialTime,0},0},{{finalTime,10},1,100}, ControlPlacement->Top]

Implementation and simulations 
Time varying parameter model

Manipulate[Module[{condizioniOK$},condizioniOK$=r1>r2;
If[condizioniOK$,Plot[{(b Linf r2 (r2 Cos[ω *t]+ω Sin[ω *t]) Exp[r2 *t]-b Linf r2^2+(L0 r2-Linf r1+Linf r1 Exp[r2 *t]) (r2^2+ω^2))*Exp[-r2* t]/(r2 (r2^2+ω^2))
},{t,initialTime,finalTime},PlotLabels->"L(t)",PlotRange->{0,All},AxesLabel->TraditionalForm/@{t,"L(t)"},LabelStyle->{FontFamily->"Times",FontSize->14},AxesStyle->Arrowheads[{0.,0.03}]],Style[Row[{"⚠️ Conditions not respected: \!\(\*SubscriptBox[\(r\), \(1\)]\)  >  \!\(\*SubscriptBox[\(r\), \(2\)]\) > 0 "," \!\(\*FractionBox[SubscriptBox[\(r\), \(1\)], SubscriptBox[\(r\), \(2\)]]\) = ",r1/r2}],Red,14,Bold]]],{{r1,0.5501},0.0001,1,0.01},{{r2,0.0601},0.0001,1,0.01},{{Linf,5.7601},0.0001,10,0.01},{{L0,3.08},0.1,10,0.01},{{initialTime,0},0},{{ω,0.1},0.0001,1,0.01},{{b,0.1},0.0001,1,0.01}, {{finalTime,100.},1,100},Initialization:>{initialTime=0,finalTime=90}, ControlPlacement->Top]

Implementation and simulations 
Time varying parameter model

KMIN[t_,L0_,L∞_, r1_,r2_, α_,w_,b_,δ_,k0_,s_,A_]:=((((-1))^(1/(1-α))) (E^((t (r2 α+δ))/(-1+α))) (k0^(α/(-1+α))) (L0^(α/(-1+α))) (r2^((2 α)/(-1+α))) ((((r2^2)+(w^2)))^(α/(-1+α))) (((b L∞ (r2^2)+(L∞ r1-L0 r2) ((r2^2)+(w^2))))^(1/(-1+α))) (((L0 ((r2^2)+(w^2))-b L∞ (r2+(Sqrt[(r2^2)+(w^2)]))))^(α/(-1+α))) (((((-1+(E^(r2 t))) L∞ r1+L0 r2) ((r2^2)+(w^2))-b L∞ r2 (r2+(E^(r2 t)) (Sqrt[(r2^2)+(w^2)]))))^(α/(-1+α))) (((r2-δ))^(1/(-1+α))) ((((r2^α) (-A (E^(t (r2 α+δ))) (k0^α) (L0^α) (r2^α) s ((((r2^2)+(w^2)))^α) (((L0 ((r2^2)+(w^2))-b L∞ (r2+(Sqrt[(r2^2)+(w^2)]))))^α) (((((-1+(E^(r2 t))) L∞ r1+L0 r2))^2) ((((r2^2)+(w^2)))^2)-2 b L∞ r2 ((-1+(E^(r2 t))) L∞ r1+L0 r2) ((r2^2)+(w^2)) (r2+(E^(r2 t)) (Sqrt[(r2^2)+(w^2)]))+(b^2) ((L∞)^2) (r2^2) ((1+(E^(2 r2 t))) (r2^2)+(E^(2 r2 t)) (w^2)+2 (E^(r2 t)) r2 (Sqrt[(r2^2)+(w^2)]))) Hypergeometric2F1[1,((r2+δ-α δ)/r2),((r2 α+δ-α δ)/r2),(((E^(r2 t)) L∞ (-b r2 (Sqrt[(r2^2)+(w^2)])+r1 ((r2^2)+(w^2))))/(b L∞ (r2^2)+(L∞ r1-L0 r2) ((r2^2)+(w^2))))]+(E^(t (r2+α δ))) r2 (((((-1+(E^(r2 t))) L∞ r1+L0 r2) ((r2^2)+(w^2))-b L∞ r2 (r2+(E^(r2 t)) (Sqrt[(r2^2)+(w^2)]))))^α) (k0 L0 ((r2^2)+(w^2)) (-b L∞ (r2^2)-(L∞ r1-L0 r2) ((r2^2)+(w^2))) (((L0 ((r2^2)+(w^2))-b L∞ (r2+(Sqrt[(r2^2)+(w^2)]))))^α) (r2-δ)+A (k0^α) (L0^α) r2 s ((((r2^2)+(w^2)))^α) ((L0^2) ((((r2^2)+(w^2)))^2)-2 b L0 L∞ ((r2^2)+(w^2)) (r2+(Sqrt[(r2^2)+(w^2)]))+(b^2) ((L∞)^2) (2 (r2^2)+(w^2)+2 r2 (Sqrt[(r2^2)+(w^2)]))) Hypergeometric2F1[1,((r2+δ-α δ)/r2),((r2 α+δ-α δ)/r2),-((b L∞ r2 (Sqrt[(r2^2)+(w^2)])-L∞ r1 ((r2^2)+(w^2)))/(b L∞ (r2^2)+(L∞ r1-L0 r2) ((r2^2)+(w^2))))]))))^(1/(1-α))))/(-b L∞ (r2^2)-L∞ r1 (r2^2)+(E^(r2 t)) L∞ r1 (r2^2)+L0 (r2^3)-L∞ r1 (w^2)+(E^(r2 t)) L∞ r1 (w^2)+L0 r2 (w^2)+b (E^(r2 t)) L∞ (r2^2) Cos[t w]+b (E^(r2 t)) L∞ r2 w Sin[t w]); KMIN[t,L0,L∞, r1,r2,α,w,b,δ,k0,s,A];
KMAX[t_,L0_,L∞_, r1_,r2_, α_,w_,b_,δ_,k0_,s_,A_]:= ((((-1))^(1/(1-α))) (E^((t (r2 α+δ))/(-1+α))) (k0^(α/(-1+α))) (L0^(α/(-1+α))) (r2^(α/(-1+α))) ((((r2^2)+(w^2)))^(α/(-1+α))) (((b L∞ (r2^2)+(L∞ r1-L0 r2) ((r2^2)+(w^2))))^(1/(-1+α))) (((L0 ((r2^2)+(w^2))+b L∞ (-r2+(Sqrt[(r2^2)+(w^2)]))))^(α/(-1+α))) (((((-1+(E^(r2 t))) L∞ r1+L0 r2) ((r2^2)+(w^2))+b L∞ r2 (-r2+(E^(r2 t)) (Sqrt[(r2^2)+(w^2)]))))^(α/(-1+α))) (((r2-δ))^(1/(-1+α))) (((A (E^(t (r2+α δ))) (k0^α) (L0^α) (r2^2) s ((((r2^2)+(w^2)))^α) (((((-1+(E^(r2 t))) L∞ r1+L0 r2) ((r2^2)+(w^2))+b L∞ r2 (-r2+(E^(r2 t)) (Sqrt[(r2^2)+(w^2)]))))^α) ((L0^2) ((((r2^2)+(w^2)))^2)+2 b L0 L∞ ((r2^2)+(w^2)) (-r2+(Sqrt[(r2^2)+(w^2)]))+(b^2) ((L∞)^2) (2 (r2^2)+(w^2)-2 r2 (Sqrt[(r2^2)+(w^2)]))) Hypergeometric2F1[1,((r2+δ-α δ)/r2),((r2 α+δ-α δ)/r2),((L∞ (b r2 (Sqrt[(r2^2)+(w^2)])+r1 ((r2^2)+(w^2))))/(b L∞ (r2^2)+(L∞ r1-L0 r2) ((r2^2)+(w^2))))]+(((L0 ((r2^2)+(w^2))+b L∞ (-r2+(Sqrt[(r2^2)+(w^2)]))))^α) ((E^(t (r2+α δ))) k0 L0 r2 ((r2^2)+(w^2)) (-b L∞ (r2^2)-(L∞ r1-L0 r2) ((r2^2)+(w^2))) (((((-1+(E^(r2 t))) L∞ r1+L0 r2) ((r2^2)+(w^2))+b L∞ r2 (-r2+(E^(r2 t)) (Sqrt[(r2^2)+(w^2)]))))^α) (r2-δ)-A (E^(t (r2 α+δ))) (k0^α) (L0^α) (r2^α) s ((((r2^2)+(w^2)))^α) (((((-1+(E^(r2 t))) L∞ r1+L0 r2))^2) ((((r2^2)+(w^2)))^2)+2 b L∞ r2 ((-1+(E^(r2 t))) L∞ r1+L0 r2) ((r2^2)+(w^2)) (-r2+(E^(r2 t)) (Sqrt[(r2^2)+(w^2)]))+(b^2) ((L∞)^2) (r2^2) ((1+(E^(2 r2 t))) (r2^2)+(E^(2 r2 t)) (w^2)-2 (E^(r2 t)) r2 (Sqrt[(r2^2)+(w^2)]))) Hypergeometric2F1[1,((r2+δ-α δ)/r2),((r2 α+δ-α δ)/r2),(((E^(r2 t)) L∞ (b r2 (Sqrt[(r2^2)+(w^2)])+r1 ((r2^2)+(w^2))))/(b L∞ (r2^2)+(L∞ r1-L0 r2) ((r2^2)+(w^2))))])))^(1/(1-α))))/(-b L∞ (r2^2)-L∞ r1 (r2^2)+(E^(r2 t)) L∞ r1 (r2^2)+L0 (r2^3)-L∞ r1 (w^2)+(E^(r2 t)) L∞ r1 (w^2)+L0 r2 (w^2)+b (E^(r2 t)) L∞ (r2^2) Cos[t w]+b (E^(r2 t)) L∞ r2 w Sin[t w]); KMAX[t,L0,L∞, r1,r2,α,w,b,δ,k0,s,A];
grules={
r2 -> 1/100,
r1 -> 1/20,
b -> 1/50,
ω -> 1/2,
L0 -> 1/2,
L∞ -> 1,
α -> 1/3,
δ -> 1/5,
k0 -> 1, 
A -> 1,
s -> 1/4};

arules={
r1 -> 1/100,
r2-> 1/20,
b -> 1/50,
ω -> 1/2,
L0 -> 1/2,
L∞ -> 1,
α -> 1/3,
δ -> 1/5,
k0 -> 1, 
A -> 1,
s -> 1/4};
G[t_,L0_,L∞_,r1_,r2_,α_,w_,b_]:=b  L∞  r2  (r2  Cos[w  t]+w  Sin[ w t])Exp[r2  t]-b  L∞  (r2^2) +(L0 r2 - L∞  r1 +  L∞  r1   Exp[r2  t])((r2^2)+(w^2));
enne[t_,L0_,L∞_,r1_,r2_,α_,w_,b_]:= r2(b  L∞  w  (w Cos[w  t]-r2  Sin[ w t])  Exp[r2  t]  +b L∞  (r2^2)-(L0 r2 -L∞ r1 ) ((r2^2)+(w^2))   )/ G[t,L0,L∞,r1,r2,α,w,b];
Clear[ode,odeG];
ode[ t_,L0_,L∞_, r1_,r2_, α_,w_,b_,δ_,k0_,s_,A_]:= KK'[t]+(enne[t,L0,L∞,r1,r2,α,w,b]+δ)*KK[t]- s A (KK[t])^α;
ode[t,L0,L∞, r1,r2,α,ω,b,δ,k,s,A];
odeG=ode[t,L0,L∞, r1,r2,α,ω,b,δ,k,s,A]/.grules//N//Chop;
nsolG=NDSolve[{odeG==0,KK[0]==1}, KK[t],{t,0,40}];
odeA=ode[t,L0,L∞, r1,r2,α,ω,b,δ,k,s,A]/.arules//N//Chop;
nsolA=NDSolve[{odeA==0,KK[0]==1}, KK[t],{t,0,40}];

gplot=Plot[{
KMAX[t,L0,L∞, r1,r2,α,ω,b,δ,k0,s,A]/.grules//N//Chop,
KMIN[t,L0,L∞, r1,r2,α,ω,b,δ,k0,s,A]/.grules//N//Chop,
KK[t]/.nsolG
},
{t,0,40},
AxesOrigin->{0,0},
PlotRange->{{-1,45},All},
AxesLabel->{Style[t,Bold],Style["k(t)",Bold]},
AxesStyle->Arrowheads[0.03],
Ticks->{All,{1,2}}]
aplot=Plot[{
KMAX[t,L0,L∞, r1,r2,α,ω,b,δ,k0,s,A]/.arules//N//Chop,
KMIN[t,L0,L∞, r1,r2,α,ω,b,δ,k0,s,A]/.arules//N//Chop,
KK[t]/.nsolA
},
{t,0,40},
AxesOrigin->{0,0},
PlotRange->{{-1,45},All},
AxesLabel->{Style[t,Bold],Style["k(t)",Bold]},
AxesStyle->Arrowheads[0.03],
Ticks->{All,{1,2}}]

Main results

•  Constant parameter case: exact solution via Hypergeometric function (Subscript[, 2])(Subscript[F, 1]) 
•  Time - varying case: no direct closed-form solution
•  Analytical upper and lower bounds for the solution of the seasonal case
•  Numerical solution, via the use of  NDSolve[ ], lies between the upper and lower bound

Conclusions

•  Extended von Bertalanffy: flexible alternative to Malthus models
•  Closed-form solutions in the constant parameters model and analytical bounds for the solution in time varying parameter model
•  Mathematica: symbolic integration, numerical solution via NDSolve[ ], and visualization tools
•  Future work: empirical applications and multi-sectoral models

Essential Bibliography

[1] Brida, J.; Maldonado, E.  Closed form solutions to a generalization of the Solow growth model. Applied Mathematical Sciences 2007,1, 1991
[2] Mingari-Scarpello, G.; Ritelli, D.,  The Solow model improved through the logistic manpower growth law, Annali Univ. Ferrara, 2003, 49, 73–83.–2000. 
[3] Mingari-Scarpello, G.; Ritelli, D, Closed form solution of a periodically forced logistic model,  Annali U. Ferrara 2008, 54, 85–94.
[4] von Bertalanffy, L., A quantitative theory of organic growth(Inquiries on growth laws. II), Human Biology 1938, 10, 181–213.
[5] Cloern,J.E.; Nichols, F.H., A von Bertalanffy growth model with a seasonally varying coefficient, J.Fish Res. Board Can. 1978,35, 1479–1482. 
[6] Sofroniou, M.; Spaletta, G., On the construction of a new generalization of Runge–Kutta methods, Electronic Notesin Theoretical Computer Science2003, 74, 189–206. 
[7] Solow,R., A Contribution to the Theory of Economic Growth. The Quarterly Journal of Economics 1956, 70, 65–94.
[8] Swan, T., Economic Growth and Capital Accumulation, Economic Record 1956, 32, 334–361.
[9] Bargellini, Ritelli, Spaletta., Extended von Bertalanffy Equation in Solow Growth Modelling. Algorithms 18, 565, 1–15, 2025
[10] Donghan, C., An improved Solow-Swan model. Chinese Quarterly Journal of Mathematics, 1998, 13, pp. 72--78
[11] Essington, T. E., Kitchell, J. F., & Walters, C. J. (2011). The von Bertalanffy growth function, bioenergetics, and the consumption rates of fish. Canadian Journal of Fisheries and Aquatic Sciences, 58, 2129–2138.
[12] Katsanevakis, S., & Maravelias, C. D. (2008). Modelling fish growth: multi-model inference as a better alternative to a priori using von Bertalanffy equation. Fish and Fisheries, 9, 178–187.
[13] Somers, I. F. (1988). On a seasonally oscillating growth functionTechnical Report, Division of Fisheries Research, CSIRO Marine Laboratories, Cleveland, Australia.

Essential Bio

Dr.Bargellini holds an Economics & Business Bachelor (Florence University 2020) and a Quantitative Finance M.Sc. (Bologna Univ. 2023). He won an Erasmus scholarship (Athens Univ. of Economics & Business 2019-2020) and a PhD scholarship (Bologna Statistics Department 2023-current). He successfully completed an Intensive Programme on Advances Topics in Energy Markets & Finance, an Advanced Training Course in Actuarial Sciences. He worked in private firms as Junior Analyst and Accounting Assistant.