--- trunk/dyn3d/coefpoly.f90 2014/03/05 14:38:41 81 +++ trunk/dyn3d/coefpoly.f 2015/02/05 12:41:08 123 @@ -1,39 +1,50 @@ +module coefpoly_m -! $Header: /home/cvsroot/LMDZ4/libf/dyn3d/coefpoly.F,v 1.1.1.1 2004/05/19 -! 12:53:05 lmdzadmin Exp $ - -SUBROUTINE coefpoly(xf1, xf2, xprim1, xprim2, xtild1, xtild2, a0, a1, a2, a3) IMPLICIT NONE - ! ... Auteur : P. Le Van ... +contains + + SUBROUTINE coefpoly(xf1, xf2, xprim1, xprim2, xtild1, xtild2, a0, a1, a2, a3) + + ! From LMDZ4/libf/dyn3d/coefpoly.F, version 1.1.1.1 2004/05/19 12:53:05 + + ! Author: P. Le Van + + ! Calcul des coefficients a0, a1, a2, a3 du polynôme de degré 3 qui + ! satisfait aux 4 équations suivantes : + + ! a0 + a1 * xtild1 + a2 * xtild1**2 + a3 * xtild1**3 = Xf1 + ! a0 + a1 * xtild2 + a2 * xtild2**2 + a3 * xtild2**3 = Xf2 + ! a1 + 2. * a2 * xtild1 + 3. * a3 * xtild1**2 = Xprim1 + ! a1 + 2. * a2 * xtild2 + 3. * a3 * xtild2**2 = Xprim2 + + ! (passe par les points (Xf(it), xtild(it)) et (Xf(it + 1), + ! xtild(it + 1)) + ! On en revient à resoudre un système de 4 équations à 4 inconnues + ! a0, a1, a2, a3. - ! Calcul des coefficients a0, a1, a2, a3 du polynome de degre 3 qui - ! satisfait aux 4 equations suivantes : + DOUBLE PRECISION, intent(in):: xf1, xf2, xprim1, xprim2, xtild1, xtild2 + DOUBLE PRECISION, intent(out):: a0, a1, a2, a3 - ! a0 + a1*xtild1 + a2*xtild1*xtild1 + a3*xtild1*xtild1*xtild1 = Xf1 - ! a0 + a1*xtild2 + a2*xtild2*xtild2 + a3*xtild2*xtild2*xtild2 = Xf2 - ! a1 + 2.*a2*xtild1 + 3.*a3*xtild1*xtild1 = Xprim1 - ! a1 + 2.*a2*xtild2 + 3.*a3*xtild2*xtild2 = Xprim2 + ! Local: + DOUBLE PRECISION xtil1car, xtil2car, derr, x1x2car - ! On en revient a resoudre un systeme de 4 equat.a 4 inconnues a0,a1,a2,a3 + !------------------------------------------------------------ - DOUBLE PRECISION xf1, xf2, xprim1, xprim2, xtild1, xtild2, xi - DOUBLE PRECISION xfout, xprim - DOUBLE PRECISION a1, a2, a3, a0, xtil1car, xtil2car, derr, x1x2car + xtil1car = xtild1 * xtild1 + xtil2car = xtild2 * xtild2 - xtil1car = xtild1*xtild1 - xtil2car = xtild2*xtild2 + derr = 2. * (xf2-xf1)/(xtild1-xtild2) - derr = 2.*(xf2-xf1)/(xtild1-xtild2) + x1x2car = (xtild1-xtild2) * (xtild1-xtild2) - x1x2car = (xtild1-xtild2)*(xtild1-xtild2) + a3 = (derr+xprim1+xprim2)/x1x2car + a2 = (xprim1-xprim2+3. * a3 * (xtil2car-xtil1car))/(2. * (xtild1-xtild2)) - a3 = (derr+xprim1+xprim2)/x1x2car - a2 = (xprim1-xprim2+3.*a3*(xtil2car-xtil1car))/(2.*(xtild1-xtild2)) + a1 = xprim1 - 3. * a3 * xtil1car - 2. * a2 * xtild1 + a0 = xf1 - a3 * xtild1 * xtil1car - a2 * xtil1car - a1 * xtild1 - a1 = xprim1 - 3.*a3*xtil1car - 2.*a2*xtild1 - a0 = xf1 - a3*xtild1*xtil1car - a2*xtil1car - a1*xtild1 + END SUBROUTINE coefpoly - RETURN -END SUBROUTINE coefpoly +end module coefpoly_m