subrutinas solve y factor

Dudas y comentarios sobre otros lenguajes de programación. Si algún lenguaje recibe suficientes preguntas le añadimos nueva categoría.
Responder
Mensaje
Autor
triff
Mensajes: 1
Registrado: 13/03/2010 2:49 pm

subrutinas solve y factor

#1 Mensaje por triff » 13/03/2010 2:57 pm

Hola a todos!! Estoy empezando a trabajar en c y tengo que
traducir de fortran a C los siguientes codigos. Para alguien
que controle un poco seguro que es algo muy facil y me ahorraria mucho tiempo. Y asi quedan escritos para el uso y disfrute de mas gente del foro, que veo que estan sedientos de metodos numericos. Estas subrutinas permiten resolver sistemas de ecuaciones lineales, o el
tipico "resuelve Ax=b".
El codigo es el siguiente:

subroutine factor(a,n)
!donde a es una matriz nxn

real(8)::a(n,n),am
integer::n,i,j,k
!calcula la descomposición LU de la matriz A
!el algoritmo de Crout
!avanza por la diagonal
do j=1,n
!calcula l(j,j) a l(n,j)
do i=j,n
do k=1,j-1
a(i,j)=a(i,j)-a(i,k)*a(k,j)
end do
end do

!calcula u(j,j+1) a u (j,n)
am=1d0/a(j,j)
do i=j+1,n
do k=1,j-1
a(j,i)=a(j,i)-a(j,k)*a(k,i)
end do
a(j,i)=a(j,i)*am
end do
end do

return
end subroutine

!esta surutina permite descomponer cualquier matriz A !en dos matrices diagonales LU y almacenarlas en la !misma matriz A. Si partimos de un sistema A*x=b, !obtendriamos LUx=b. Resolvemos con la siguiente !subrutina haciendo Ly=b , y Ux=y.

subroutine solve(a,b,n)
real(8)::a(n,n),b(n)
integer::n,i,j
!resuelve el sistema Ax=b con la matriz A factorizada en !LU mediante el algoritmo de Crout

!resuelve Ly=b
do i=1,n
do j=1,i-1
b(i)=b(i)-a(i,j)*b(j)
end do
b(i)=b(i)/a(i,i)
end do

!resuelve Ux=y
do i=n,1,-1
do j=i+1,n
b(i)=b(i)-a(i,j)*b(i)
end do
end do

return
end subroutine




he probado a intentar pasar una de las subrutinas a c, pero con
mis conocimientos tan escasos de c...

void factor(int a[n][n],int n) {
/* le estamos poniendo int para provar */
int am,i,j,k;

/*calcula la descomposicion LU de la matriz A*/
/*el algoritmo de Crout*/
/* avanza por la diagonal*/
/*calcula l(i,j) a l(n,j)*/
for (j=1;j<=n;j++)
{
for (i=j;i<=n;i++)
{
for (k=1;k<=j-1;k++)
{
a(i,j)=a(i,j)-a(i,k)*a(k,j);
}
}
}
/*calcula u(j,j+1) a u(j,n) */
am=1/a(j,j);
for (i=j+1;i<=n;i++)
{
for (k=1;k<=j-1;k++)
{
a(j,i)=a(j,i)-a(j,k)*a(k,i);
}
}
a(j,i)=a(j,i)*am;
}
}

Responder

¿Quién está conectado?

Usuarios navegando por este Foro: No hay usuarios registrados visitando el Foro y 1 invitado