inversion_solver Subroutine

public subroutine inversion_solver(n, rhs, x, tmat)

Uses

  • proc~~inversion_solver~~UsesGraph proc~inversion_solver inversion_solver module~mod_memory mod_memory proc~inversion_solver->module~mod_memory module~mod_io mod_io module~mod_memory->module~mod_io iso_c_binding iso_c_binding module~mod_memory->iso_c_binding module~mod_constants mod_constants module~mod_memory->module~mod_constants module~mod_io->module~mod_constants module~mod_constants->iso_c_binding

Solve the linear system directly inverting the matrix: This is highly unefficient and should only be used for testing other methods of solution.

Arguments

Type IntentOptional Attributes Name
integer(kind=ip), intent(in) :: n

Size of the matrix
real(kind=rp), intent(in), dimension(n) :: rhs

Right hand side of the linear system
real(kind=rp), intent(out), dimension(n) :: x

In output the solution of the linear system
real(kind=rp), intent(in), dimension(n, n) :: tmat

Polarization matrix TODO

Calls

proc~~inversion_solver~~CallsGraph proc~inversion_solver inversion_solver interface~mallocate mallocate proc~inversion_solver->interface~mallocate dgetrf dgetrf proc~inversion_solver->dgetrf dgetri dgetri proc~inversion_solver->dgetri dgemm dgemm proc~inversion_solver->dgemm interface~mfree mfree proc~inversion_solver->interface~mfree proc~r_alloc1 r_alloc1 interface~mallocate->proc~r_alloc1 proc~i_alloc2 i_alloc2 interface~mallocate->proc~i_alloc2 proc~i_alloc1 i_alloc1 interface~mallocate->proc~i_alloc1 proc~r_alloc3 r_alloc3 interface~mallocate->proc~r_alloc3 proc~r_alloc2 r_alloc2 interface~mallocate->proc~r_alloc2 proc~i_alloc3 i_alloc3 interface~mallocate->proc~i_alloc3 proc~l_alloc1 l_alloc1 interface~mallocate->proc~l_alloc1 proc~l_alloc2 l_alloc2 interface~mallocate->proc~l_alloc2 proc~l_free2 l_free2 interface~mfree->proc~l_free2 proc~l_free1 l_free1 interface~mfree->proc~l_free1 proc~r_free1 r_free1 interface~mfree->proc~r_free1 proc~i_free1 i_free1 interface~mfree->proc~i_free1 proc~i_free3 i_free3 interface~mfree->proc~i_free3 proc~i_free2 i_free2 interface~mfree->proc~i_free2 proc~r_free3 r_free3 interface~mfree->proc~r_free3 proc~r_free2 r_free2 interface~mfree->proc~r_free2 proc~chk_free chk_free proc~l_free2->proc~chk_free proc~l_free1->proc~chk_free proc~chk_alloc chk_alloc proc~r_alloc1->proc~chk_alloc proc~memory_init memory_init proc~r_alloc1->proc~memory_init proc~i_alloc2->proc~chk_alloc proc~i_alloc2->proc~memory_init proc~i_alloc1->proc~chk_alloc proc~i_alloc1->proc~memory_init proc~r_alloc3->proc~chk_alloc proc~r_alloc3->proc~memory_init proc~r_alloc2->proc~chk_alloc proc~r_alloc2->proc~memory_init proc~i_alloc3->proc~chk_alloc proc~i_alloc3->proc~memory_init proc~l_alloc1->proc~chk_alloc proc~l_alloc1->proc~memory_init proc~l_alloc2->proc~chk_alloc proc~l_alloc2->proc~memory_init proc~r_free1->proc~chk_free proc~i_free1->proc~chk_free proc~i_free3->proc~chk_free proc~i_free2->proc~chk_free proc~r_free3->proc~chk_free proc~r_free2->proc~chk_free proc~fatal_error fatal_error proc~chk_free->proc~fatal_error proc~chk_alloc->proc~fatal_error proc~ommp_message ommp_message proc~fatal_error->proc~ommp_message proc~close_output close_output proc~fatal_error->proc~close_output proc~close_output->proc~ommp_message

Called by

proc~~inversion_solver~~CalledByGraph proc~inversion_solver inversion_solver proc~polarization polarization proc~polarization->proc~inversion_solver proc~ommp_get_polelec_energy ommp_get_polelec_energy proc~ommp_get_polelec_energy->proc~polarization proc~ommp_set_external_field ommp_set_external_field proc~ommp_set_external_field->proc~polarization proc~polelec_geomgrad polelec_geomgrad proc~polelec_geomgrad->proc~polarization proc~ommp_get_full_ele_energy ommp_get_full_ele_energy proc~ommp_get_full_ele_energy->proc~ommp_get_polelec_energy program~test_si_potential test_SI_potential program~test_si_potential->proc~ommp_get_polelec_energy program~test_si_potential->proc~ommp_set_external_field proc~c_ommp_set_external_field C_ommp_set_external_field proc~c_ommp_set_external_field->proc~ommp_set_external_field proc~ommp_set_external_field_nomm ommp_set_external_field_nomm proc~ommp_set_external_field_nomm->proc~ommp_set_external_field proc~c_ommp_get_polelec_energy C_ommp_get_polelec_energy proc~c_ommp_get_polelec_energy->proc~ommp_get_polelec_energy proc~c_ommp_set_external_field_nomm C_ommp_set_external_field_nomm proc~c_ommp_set_external_field_nomm->proc~ommp_set_external_field proc~ommp_polelec_geomgrad ommp_polelec_geomgrad proc~ommp_polelec_geomgrad->proc~polelec_geomgrad proc~ommp_full_geomgrad ommp_full_geomgrad proc~ommp_full_geomgrad->proc~polelec_geomgrad proc~ommp_get_full_energy ommp_get_full_energy proc~ommp_get_full_energy->proc~ommp_get_full_ele_energy proc~c_ommp_polelec_geomgrad C_ommp_polelec_geomgrad proc~c_ommp_polelec_geomgrad->proc~ommp_polelec_geomgrad proc~ommptest_fakeqm_internal_geomgrad ommptest_fakeqm_internal_geomgrad proc~ommptest_fakeqm_internal_geomgrad->proc~ommp_full_geomgrad proc~c_ommp_full_geomgrad C_ommp_full_geomgrad proc~c_ommp_full_geomgrad->proc~ommp_full_geomgrad proc~c_ommp_get_full_ele_energy C_ommp_get_full_ele_energy proc~c_ommp_get_full_ele_energy->proc~ommp_get_full_ele_energy proc~ommptest_totalqmmm_geomgrad ommptest_totalqmmm_geomgrad proc~ommptest_totalqmmm_geomgrad->proc~ommp_full_geomgrad proc~ommptest_fakeqm_linkatom_geomgrad ommptest_fakeqm_linkatom_geomgrad proc~ommptest_fakeqm_linkatom_geomgrad->proc~ommp_full_geomgrad proc~c_ommp_get_full_energy C_ommp_get_full_energy proc~c_ommp_get_full_energy->proc~ommp_get_full_energy

Contents

Source Code

inversion_solver

Source Code

    subroutine inversion_solver(n, rhs, x, tmat)
        !! Solve the linear system directly inverting the matrix:
        !! $$\mathbf A \mathbf x = \mathbf B $$
        !! $$ \mathbf x = \mathbf A ^-1 \mathbf B $$
        !! This is highly unefficient and should only be used for testing 
        !! other methods of solution.

        use mod_memory, only: mallocate, mfree
        
        implicit none
        
        integer(ip), intent(in) :: n
        !! Size of the matrix
        real(rp), dimension(n), intent(in) :: rhs
        !! Right hand side of the linear system
        real(rp), dimension(n), intent(out) :: x
        !! In output the solution of the linear system
        real(rp), dimension(n, n), intent(in) :: tmat
        !! Polarization matrix TODO

        integer(ip) :: info
        integer(ip), dimension(:), allocatable :: ipiv
        real(rp), dimension(:), allocatable :: work
        real(rp), dimension(:,:), allocatable :: TMatI
        
        call mallocate('inversion_solver [TMatI]', n, n, TMatI)
        call mallocate('inversion_solver [work]', n, work)
        call mallocate('inversion_solver [ipiv]', n, ipiv)
        
        ! Initialize inverse polarization matrix
        TMatI = TMat
        
        !Compute the inverse of TMat
        call dgetrf(n, n, TMatI, n, iPiv, info)
        call dgetri(n, TMatI, n, iPiv, Work, n, info)
        
        ! Calculate dipoles with matrix inversion
        call dgemm('N', 'N', n, 1, n, 1.0_rp, TMatI, n, rhs, n, 0.0_rp, x, n)
        
        call mfree('inversion_solver [TMatI]', TMatI)
        call mfree('inversion_solver [work]', work)
        call mfree('inversion_solver [ipiv]', ipiv)
      
    end subroutine inversion_solver