! mppt_4pam.f90
! Cálculo de MPPT para un módulo FV con modelo de 4 parámetros (basado en Type103)
! - Ajusta Rs por bisección (series4pam) con datos de catálogo y coeficientes de temperatura
! - Calcula gam (gamma), IL_ref, Io_ref en condiciones de referencia
! - Escala IL e Io a (G, Tc)
! - Calcula Voc, Isc
! - Encuentra (Imp, Vmp, Pmp) con Newton (condición dP/dI = 0) como en el original
! Autor: (tu nombre)
program mppt_4pam
implicit none
!==================== Entrada del usuario ====================
double precision :: G, Tc_C, Tref_C, Itot_ref
double precision :: Isc_ref, Voc_ref, Imp_ref, Vmp_ref
double precision :: muIsc_ref, muVoc_ref, eg_eV
double precision :: sercell, area_m2
!==================== Variables internas =====================
double precision :: Tc, Tref, q_bz, gam, a, Rs, Rsh
double precision :: IL_ref, Io_ref, IL, Io, Voc, Isc
double precision :: imp, vmp, pmp, ff, eta
integer :: ierr
! ======= Parámetros de ejemplo (edítalos a tus datos) =======
Isc_ref = 8.2d0 ! [A] corriente de cortocircuito a ref
Voc_ref = 37.0d0 ! [V] voltaje de circuito abierto a ref
Imp_ref = 7.6d0 ! [A] corriente en MPP a ref
Vmp_ref = 30.0d0 ! [V] voltaje en MPP a ref
Tref_C = 25.0d0 ! [°C] temperatura de referencia (STC)
Itot_ref = 1000.d0 ! [W/m2] irradiancia de referencia
muIsc_ref = 0.004d0 ! [A/K] coef. temp. de Isc (ejemplo)
muVoc_ref = -0.12d0 ! [V/K] coef. temp. de Voc (ejemplo)
eg_eV = 1.12d0 ! [eV] energía de banda prohibida
sercell = 60.d0 ! [-] celdas en serie en el módulo
area_m2 = 1.94d0 ! [m2] área del módulo (ejemplo)
! ======= Condiciones de operación (edita a tu caso) ========
G = 850.d0 ! [W/m2] irradiancia sobre el módulo
Tc_C = 45.d0 ! [°C] temperatura de célula (si no la sabes, usa NOCT o balance aparte)
! ======= Conversión de unidades y constantes ================
Tc = Tc_C + 273.15d0
Tref = Tref_C + 273.15d0
q_bz = 11604.45d0 ! q/k en [K/V] (mismo del Type103)
Rsh = 1.d10 ! Muy grande para el 4-parámetros (sin derivación)
! ======= 1) Calcular Rs por bisección (como Type103) ========
call series4pam_103(Rs, q_bz, sercell, Tref, Imp_ref, Vmp_ref, Isc_ref, Voc_ref, muVoc_ref, muIsc_ref, eg_eV, ierr)
if (ierr /= 0) then
write(*,*) 'Aviso: Rs muy pequeña; revisa parámetros de catálogo.'
end if
! ======= 2) gamma (gam) a referencia y parámetros base ======
gam = q_bz*(Vmp_ref - Voc_ref + Imp_ref*Rs)/(Tref*dlog(1.d0 - Imp_ref/Isc_ref))
a = gam/sercell
IL_ref = Isc_ref
Io_ref = IL_ref/dexp(q_bz*Voc_ref/(gam*Tref))
! ======= 3) Escalamiento a (G, Tc) ==========================
IL = (G/Itot_ref)*(IL_ref + muIsc_ref*(Tc - Tref))
if (IL < 0.d0) IL = 0.d0
Io = Io_ref * ( (Tc/Tref)**3 ) * dexp( (q_bz*eg_eV/a)*((1.d0/Tref) - (1.d0/Tc)) )
if (IL > 0.d0) then
Voc = gam*Tc/q_bz * dlog(IL/Io + 1.d0)
else
Voc = 0.d0
end if
Isc = IL
! ======= 4) MPPT por Newton en I (dP/dI = 0) ================
call mppt_newton_4pam(IL, Io, gam, Tc, q_bz, Rs, Imp_ref, muIsc_ref, Itot_ref, G, Tref, imp, vmp, pmp)
! ======= 5) Métricas adicionales ============================
if (Voc > 0.d0 .and. Isc > 0.d0) then
ff = (vmp*imp)/(Voc*Isc)
else
ff = 0.d0
end if
if (G > 1.d-9 .and. area_m2 > 0.d0) then
eta = pmp / (G*area_m2)
else
eta = 0.d0
end if
! ======= Salidas ============================================
write(*,'(a)') '========== MPPT (modelo FV 4-parámetros) =========='
write(*,'(a,f8.2)') 'Irradiancia [W/m2] : ', G
write(*,'(a,f8.2)') 'Tc [°C] : ', Tc_C
write(*,'(a,f8.3)') 'Rs [ohm] : ', Rs
write(*,'(a,f8.3)') 'gamma [-] : ', gam
write(*,'(a,f8.3)') 'IL [A] : ', IL
write(*,'(a,1pe12.4)') 'Io [A] : ', Io
write(*,'(a,f8.3)') 'Voc [V] : ', Voc
write(*,'(a,f8.3)') 'Isc [A] : ', Isc
write(*,'(a,f8.3)') 'Vmp [V] : ', vmp
write(*,'(a,f8.3)') 'Imp [A] : ', imp
write(*,'(a,f8.3)') 'Pmp [W] : ', pmp
write(*,'(a,f8.3)') 'FF [-] : ', ff
write(*,'(a,f8.3)') 'Eta [-] : ', eta
contains
! ------------------------------------------------------------
! MPPT por Newton (idéntico en forma al Type103):
! - Variable de iteración: i (corriente en el MPP)
! - Condición: dP/di = 0 -> f1(i)=0 con derivada f1p(i)
! - Tras converger, vmp = ln(1 + (IL-imp)/Io)*(T*gam/q_bz) - imp*Rs
! ------------------------------------------------------------
subroutine mppt_newton_4pam(IL, Io, gam, T, q_bz, Rs, Imp_ref, muIsc_ref, Itot_ref, G, Tref, imp, vmp, pmp)
implicit none
double precision, intent(in) :: IL, Io, gam, T, q_bz, Rs
double precision, intent(in) :: Imp_ref, muIsc_ref, Itot_ref, G, Tref
double precision, intent(out) :: imp, vmp, pmp
double precision :: imxo, imxn, f1, f1p
integer :: it
! Arranque: corriente MPP escalada por irradiancia y temperatura (como en Type103)
if (G > 0.d0) then
imxn = (G/Itot_ref)*(Imp_ref + muIsc_ref*(T - Tref))
else
imp = 0.d0; vmp = 0.d0; pmp = 0.d0
return
end if
do it = 1, 200
imxo = imxn
! No permitir i > IL (corriente no puede exceder a la fotogenerada)
if (imxo > IL) imxo = 0.99d0*IL
if (Io > 0.d0) then
f1 = imxo + (imxo - IL - Io) * ( dlog((IL - imxo + Io)/Io) - imxo*q_bz*Rs/(gam*T) ) &
/ ( 1.d0 + (IL - imxo + Io)*(q_bz*Rs/(gam*T)) )
f1p = 2.d0 + ( dlog((IL - imxo + Io)/Io) - (q_bz*Rs*imxo/(gam*T)) ) &
/ ( (1.d0 + (IL - imxo + Io)*(q_bz*Rs/(gam*T)))**2 )
else
f1 = 0.d0
f1p = 1.d0
end if
imxn = imxo - f1/f1p
if (dabs(imxn - imxo) <= 1.d-3) exit
end do
imp = max(0.d0, min(imxn, IL))
if (Io > 0.d0) then
vmp = dlog(1.d0 + (IL - imp)/Io)*(T*gam/q_bz) - imp*Rs
else
vmp = 0.d0
end if
pmp = imp*vmp
end subroutine mppt_newton_4pam
! ------------------------------------------------------------
! Cálculo de Rs (bisección) siguiendo el Type103.
! Usa consistencia con coeficientes de temperatura para Voc e Isc.
! ------------------------------------------------------------
subroutine series4pam_103(RSeries,q_bz,sercell,Tcell_ref,Imppt_ref,Vmppt_ref,Isc_ref,Voc_ref,muVoc_ref,muIsc_ref,eg,ierr)
implicit none
double precision, intent(out) :: RSeries
double precision, intent(in) :: q_bz, sercell, Tcell_ref
double precision, intent(in) :: Imppt_ref, Vmppt_ref, Isc_ref, Voc_ref
double precision, intent(in) :: muVoc_ref, muIsc_ref, eg
double precision :: aup, alow, anew, gamup, gamlow, gamnew
double precision :: rsup, rslow, rsnew, ioup, iolow, ionew
double precision :: fup, flow, fnew
integer :: it
integer :: ierr
ierr = 0
! Límite superior (derivado del Type103)
rsup = ((sercell*Tcell_ref*dlog(1.d0 - Imppt_ref/Isc_ref)/q_bz)+Voc_ref - Vmppt_ref)/Imppt_ref
aup = 1.d0
gamup = sercell
ioup = Isc_ref*dexp(-q_bz*Voc_ref/(gamup*Tcell_ref))
! Límite inferior
rslow = 0.d0
gamlow = q_bz*(Vmppt_ref - Voc_ref)/(Tcell_ref*dlog(1.d0 - Imppt_ref/Isc_ref))
alow = gamlow/sercell
iolow = Isc_ref*dexp(-q_bz*Voc_ref/(gamlow*Tcell_ref))
do it = 1, 10000
if (dabs(rsup - rslow) <= 5.d-4) exit
rsnew = 0.5d0*(rsup + rslow)
gamnew = q_bz*(Vmppt_ref - Voc_ref + Imppt_ref*rsnew)/(Tcell_ref*dlog(1.d0 - Imppt_ref/Isc_ref))
anew = gamnew/sercell
ionew = Isc_ref*dexp(-q_bz*Voc_ref/(gamnew*Tcell_ref))
fup = -muVoc_ref + (gamup/q_bz)*( dlog(1.d0 + Isc_ref/ioup) + (Tcell_ref/(Isc_ref + ioup)) &
*(muIsc_ref - Isc_ref*((q_bz*eg/(aup*Tcell_ref*Tcell_ref)) + 3.d0/Tcell_ref)) )
flow = -muVoc_ref + (gamlow/q_bz)*( dlog(1.d0 + Isc_ref/iolow) + (Tcell_ref/(Isc_ref + iolow)) &
*(muIsc_ref - Isc_ref*((q_bz*eg/(alow*Tcell_ref*Tcell_ref)) + 3.d0/Tcell_ref)) )
fnew = -muVoc_ref + (gamnew/q_bz)*( dlog(1.d0 + Isc_ref/ionew) + (Tcell_ref/(Isc_ref + ionew)) &
*(muIsc_ref - Isc_ref*((q_bz*eg/(anew*Tcell_ref*Tcell_ref)) + 3.d0/Tcell_ref)) )
if (flow*fnew < 0.d0) then
rsup = rsnew
gamup = q_bz*(Vmppt_ref - Voc_ref + Imppt_ref*rsup)/(Tcell_ref*dlog(1.d0 - Imppt_ref/Isc_ref))
aup = gamup/sercell
ioup = Isc_ref*dexp(-q_bz*Voc_ref/(gamup*Tcell_ref))
else
rslow = rsnew
gamlow = q_bz*(Vmppt_ref - Voc_ref + Imppt_ref*rslow)/(Tcell_ref*dlog(1.d0 - Imppt_ref/Isc_ref))
alow = gamlow/sercell
iolow = Isc_ref*dexp(-q_bz*Voc_ref/(gamlow*Tcell_ref))
end if
end do
RSeries = rsnew
if (RSeries < 1.d-3) ierr = 1
end subroutine series4pam_103
end program mppt_4pam
To embed this program on your website, copy the following code and paste it into your website's HTML: