Por que meu solucionador linear iterativo não está convergindo?

26

O que pode dar errado ao usar métodos Krylov pré-condicionados da KSP ( pacote de solucionadores lineares do PETSc ) para resolver um sistema linear esparso, como os obtidos pela discretização e linearização de equações diferenciais parciais?

Que etapas posso tomar para determinar o que está errado com o meu problema?

Que mudanças posso fazer para resolver com êxito e eficiência meu sistema linear?

Jed Brown
fonte
Você pretende que essa pergunta seja uma pergunta sobre solucionadores lineares iterativos especificamente no PETSc (que é o que eu teria coletado do texto do corpo da pergunta) ou uma pergunta sobre possíveis falhas algorítmicas de solucionadores lineares iterativos em um software principalmente - contexto agnóstico (que é o que eu teria reunido ao olhar apenas para a manchete)?
Geoff Oxberry
4
Tem a petscetiqueta. A metodologia é geral, mas acho que a resposta seria menos útil se cada "tente isso" também não incluísse o "como". Como alternativa, o "como" precisaria ser muito mais longo (e mais suscetível a erros para o visualizador) se precisasse ser explicado de uma maneira independente de software. Se alguém quiser explicar como fazer todas essas coisas usando um pacote diferente, felizmente tornarei a questão independente de software e mudarei minha resposta para indicar que ela descreve o que fazer no PETSc. Nota: eu adicionei isso, que é uma versão aprimorada de uma FAQ, para que eu possa gostar de pessoas neste site.
Jed Brown

Respostas:

26

Conselho inicial

  • Sempre corra -ksp_converged_reason -ksp_monitor_true_residualao tentar descobrir por que um método não está convergindo.
  • Faça o tamanho do problema e o número de processos o mais pequeno possível para demonstrar a falha. Você geralmente obtém insight ao determinar quais pequenos problemas exibem o comportamento que está causando a quebra do método e o tempo de resposta é reduzido. Além disso, existem algumas técnicas de investigação que só podem ser usadas para sistemas pequenos.
  • Se o problema surgir apenas após um grande número de etapas, etapas de continuação ou etapas de solução não linear, considere gravar o estado do modelo quando ocorrer uma falha, para que você possa experimentar rapidamente.
  • Como alternativa, especialmente se o seu software não tiver capacidade de ponto de verificação, use -ksp_view_binaryou MatView()para salvar o sistema linear, use o código em $PETSC_DIR/src/ksp/ksp/examples/tutorials/ex10.cpara ler na matriz e resolvê-lo (possivelmente com um número diferente de processos). Isso requer uma matriz montada, portanto sua utilidade pode ser um pouco limitada.
  • Existem muitas opções possíveis de solucionador (por exemplo, um número infinito disponível na linha de comando no PETSc devido a um número arbitrário de níveis de composição); consulte esta pergunta para obter conselhos gerais sobre a escolha de solucionadores lineares.

Razões comuns para o KSP não convergir

  • As equações são singulares por acidente (por exemplo, se esqueceu de impor condições de contorno). Verifique se há um pequeno problema usando -pc_type svd -pc_svd_monitor. Tente também um solucionador direto com -pc_type lu(por meio de um pacote de terceiros em paralelo, por exemplo -pc_type lu -pc_factor_mat_solver_package superlu_dist).
  • As equações são intencionalmente singulares (por exemplo, espaço nulo constante), mas o método de Krylov não foi informado, veja KSPSetNullSpace().
  • As equações são intencionalmente singulares e KSPSetNullSpace()foram usadas, mas o lado direito não é consistente. Talvez você precise ligar MatNullSpaceRemove()do lado direito antes de ligar KSPSolve().
  • As equações são indefinidas para que os pré-condicionantes padrão não funcionem. Normalmente, você saberá disso pela física, mas pode verificar com -ksp_compute_eigenvalues -ksp_gmres_restart 1000 -pc_type none. Para problemas simples do ponto de sela, tente -pc_type fieldsplit -pc_fieldsplit_type schur -pc_fieldsplit_detect_saddle_point. Consulte o manual do usuário e a página do manual PCFIELDSPLIT para obter mais detalhes. Para problemas mais difíceis, leia a literatura para encontrar métodos robustos e pergunte aqui (ou [email protected]ou [email protected]) se deseja conselhos sobre como implementá-los. Por exemplo, veja esta pergunta para Helmholtz de alta frequência. Para tamanhos modestos de problemas, veja se você consegue conviver usando apenas um solucionador direto.
  • Se o método convergir em residual pré-condicionado, mas não em residual verdadeiro, o pré-condicionador provavelmente é singular ou quase isso. Isso é comum em problemas de ponto de sela (por exemplo, fluxo incompressível) ou operadores fortemente não simétricos (por exemplo, problemas hiperbólicos de baixo Mach com grandes intervalos de tempo).
  • O pré-condicionador é muito fraco ou é instável. Veja se -pc_type asm -sub_pc_type lumelhora a taxa de convergência. Se o GMRES estiver perdendo muito progresso na reinicialização, verifique se a ajuda é reiniciada por mais tempo -ksp_gmres_restart 300. Se uma transposição estiver disponível, tente -ksp_type bcgsou outros métodos que não exijam reinicialização. (Observe que a convergência com esses métodos é frequentemente irregular.)
  • A matriz de pré-condicionamento pode não estar próxima do operador (possivelmente não montado). Tente resolver com um solucionador direto, seja em série -pc_type luou em paralelo, usando um pacote de terceiros (por exemplo -pc_type lu -pc_factor_mat_solver_package superlu_dist, ou mumps). O método deve convergir em uma iteração se as matrizes forem iguais e em um número "pequeno" de iterações caso contrário. Tente -snes_type testverificar as matrizes se estiver resolvendo um problema não linear.
  • O pré-condicionador não é linear (por exemplo, uma solução iterativa aninhada), tente -ksp_type fgmres or -ksp_type gcr.
  • Você está usando multigrid geométrica, mas algumas equações (geralmente condições de contorno) não são dimensionadas de forma compatível entre os níveis. Tente -pc_mg_galerkinconstruir algebricamente um operador grosseiro corretamente dimensionado ou certifique-se de que todas as equações sejam dimensionadas da mesma maneira, se você desejar usar níveis grosseiros rediscretizados.
  • A matriz está muito mal condicionada. Verifique o número da condição usando os métodos descritos aqui . Tente aprimorá-lo escolhendo a escala relativa dos componentes / condições de contorno. Tente -ksp_diagonal_scale -ksp_diagonal_scale_fix. Talvez mude a formulação do problema para produzir equações algébricas mais amigáveis. Se você não conseguir corrigir a escala, pode ser necessário usar um solucionador direto.
  • A matriz é não linear (por exemplo, avaliada usando diferenciação finita de uma função não linear). Tente diferentes parâmetros de diferenciação (por exemplo -mat_mffd_type ds). Tente usar maior precisão para tornar a diferença mais precisa ./configure --with-precision=__float128 --download-f2cblaslapack,. Verifique se converge em regimes de parâmetros "mais fáceis".
  • Um método simétrico está sendo usado para um problema não simétrico.
  • O clássico Gram-Schmidt está se tornando instável, tente -ksp_gmres_modifiedgramschmidtou use um método que ortogonalize diferentemente, por exemplo -ksp_type gcr.
Jed Brown
fonte
16

Meu conselho para os alunos é tentar um solucionador direto nesses casos. A razão é que existem duas classes de razões pelas quais um solucionador pode não convergir: (i) a matriz está incorreta ou (ii) há um problema com o solucionador / pré-condicionador. Os solucionadores diretos quase sempre produzem algo que você pode comparar com a solução que você espera; portanto, se a resposta do solucionador direto parecer correta, você saberá que o problema está no solucionador / pré-condição iterativa. Por outro lado, se a resposta parecer errada, o problema está em montar a matriz e o lado direito.

Normalmente, apenas uso o UMFPACK como solucionador direto. Tenho certeza de que é simples tentar algo semelhante com o PETSC.

Wolfgang Bangerth
fonte
5
-pc_type lu -pc_factor_mat_solver_type umfpackpara usar UMFPACK (ou -pc_type cholesky -pc_factor_mat_solver_package cholmodpara problemas de SPD) através do PETSc, mas observe que UMFPACK e CHOLMOD são seriais. Para paralelo, utilização -pc_factor_mat_solver_package superlu_distou mumps, pastix, spooles.
Jed Brown
2
Só para ficar claro, a opção completa definida para o uso (por exemplo) superlu_distseria -ksp_type preonly -pc_type lu -pc_factor_mat_solver_package superlu_dist. Isso esta certo?
Leon Avery
Eu não sei. O que acontece se você fizer isso?
Wolfgang Bangerth