Início / coding / Perguntas / 79198541
Accepted
TobiR
Asked: 2024-11-18 10:21:23 +0800 CST

É seguro passar um ponteiro de uma matriz não contígua para uma rotina Lapack95 em Fortran?

  • 772

Gostaria de chamar a syevdsub-rotina Lapack95 para uma matriz não contígua da seguinte maneira:

real :: mat(15000, 15000), vec(15000)
mat=1.d0
associate(eig_vects=>mat(:10000,:10000),eig_vals=>vec(:10000))
   call syevd(eig_vects,eig_evals,'V')
end associate

É seguro fazer isso? Anteriormente, tive problemas ao passar fatias de array não contíguas para uma subrotina (o array temporário criado por essa subrotina não pôde ser alocado). Posso esperar problemas semelhantes ao usar o eig_vectsponteiro como argumento?

arrays

1 respostas

  1. Best Answer

    Lapack95 não é uma implementação independente do Lapack, mas sim um wrapper para o Lapack clássico, com uma interface mais moderna. Em particular, os argumentos de array fictícios do Lapack95 são assumidos como forma (que pode ser descontínuo), enquanto eles são assumidos como tamanho no Lapack clássico (portanto, contíguos).

    Como o fim da cadeia de chamadas (Lapack clássico) espera arrays contíguos, uma cópia temporária tem que ser feita de qualquer forma se você quiser processar arrays descontíguos. A única questão é "onde":

    • associate( eig_vects=>mat(:10000,:10000) )não faz nenhuma cópia, eigs_vectsé apenas um tipo de alias aqui
    • call syevd(eig_vects,eig_evals,'V')também não faz nenhuma cópia, pois a interface do Lapack95 syevdpode manipular uma matriz descontínua.
    • Então, a cópia temporária provavelmente acontecerá na chamada interna para a rotina Lapack *syevd; mas você não tem controle sobre como ela é realizada. Por exemplo, normalmente o compilador Intel usa a pilha para matrizes temporárias, e ele falhará com uma matriz tão grande.

    Conforme sugerido por @VladimirFГероямслава, é melhor alocar explicitamente uma matriz contígua para trabalhar de uma vez por todas:

    real, allocatable :: eigs_vect(:,:)
allocate( eigs_vect, source=mat(:10000,:10000) )

    Observe que ele vec(:10000)é contíguo, então você não precisa alocar outro array para ele.

    Note que você também pode chamar diretamente as rotinas clássicas de Lapack, já que suas interfaces podem lidar com esse tipo de situação sem nenhuma cópia, graças ao desacoplamento devido ao tamanho do problema e da "dimensão principal":

    call ssyevd( 'V', 'U', 10000, mat, 15000, vec & 
           , work, lwork, iwork, liwork, stat )
    • 1

relate perguntas