CP2K 2.4 (Revision 12889)

cp_fm_types.f90

Go to the documentation of this file.
00001 !-----------------------------------------------------------------------------!
00002 !   CP2K: A general program to perform molecular dynamics simulations         !
00003 !   Copyright (C) 2000 - 2013  CP2K developers group                          !
00004 !-----------------------------------------------------------------------------!
00005 
00006 ! *****************************************************************************
00013 MODULE cp_fm_types
00014   USE cp_array_r_utils,                ONLY: cp_2d_r_write
00015   USE cp_blacs_calls,                  ONLY: cp_blacs_dgebr2d,&
00016                                              cp_blacs_dgebs2d,&
00017                                              cp_blacs_gridexit,&
00018                                              cp_blacs_gridinfo,&
00019                                              cp_blacs_gridinit
00020   USE cp_files,                        ONLY: close_file,&
00021                                              open_file
00022   USE cp_fm_struct,                    ONLY: cp_fm_struct_equivalent,&
00023                                              cp_fm_struct_get,&
00024                                              cp_fm_struct_release,&
00025                                              cp_fm_struct_retain,&
00026                                              cp_fm_struct_type,&
00027                                              cp_fm_struct_write
00028   USE cp_para_types,                   ONLY: cp_blacs_env_type,&
00029                                              cp_para_env_type
00030   USE kinds,                           ONLY: default_path_length,&
00031                                              dp,&
00032                                              sp
00033   USE message_passing,                 ONLY: cp2k_is_parallel,&
00034                                              mp_bcast,&
00035                                              mp_max,&
00036                                              mp_recv,&
00037                                              mp_send,&
00038                                              mp_sum
00039   USE parallel_rng_types,              ONLY: UNIFORM,&
00040                                              create_rng_stream,&
00041                                              delete_rng_stream,&
00042                                              get_rng_stream,&
00043                                              random_numbers,&
00044                                              reset_to_next_rng_substream,&
00045                                              rng_stream_type
00046   USE string_utilities,                ONLY: compress
00047   USE timings,                         ONLY: timeset,&
00048                                              timestop
00049 #include "cp_common_uses.h"
00050 
00051   IMPLICIT NONE
00052 
00053   PRIVATE
00054 
00055   CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'cp_fm_types'
00056   LOGICAL, PARAMETER          :: debug_this_module = .TRUE.
00057 
00058   INTEGER, SAVE :: last_fm_id_nr = 0
00059 
00060   PUBLIC :: cp_fm_type, &
00061             cp_fm_p_type
00062 
00063   PUBLIC :: cp_fm_add_to_element,&
00064             cp_fm_create,&
00065             cp_fm_retain,&
00066             cp_fm_release,&
00067             cp_fm_get_info, &
00068             cp_fm_set_element, &
00069             cp_fm_get_element, &
00070             cp_fm_get_diag, &          ! get diagonal
00071             cp_fm_set_diag, &          ! set diagonal
00072             cp_fm_set_all, &           ! set all elements and diagonal
00073             cp_fm_set_submatrix, &     ! set a submatrix to given values
00074             cp_fm_get_submatrix, &     ! get a submatrix of given values
00075             cp_fm_init_random, &
00076             cp_fm_maxabsval, &         ! find the maximum absolute value
00077             cp_fm_maxabsrownorm, &     ! find the maximum of the sum of the abs of the elements of a row
00078             cp_fm_write, &             ! writes lots of information and data (might need further development)
00079             cp_fm_to_fm,&              ! copy (parts of) a fm to a fm
00080             cp_fm_vectorsnorm,&          ! compute the norm of the column-vectors
00081             cp_fm_to_fm_submat,&          ! copy (parts of) a fm to a fm
00082             cp_fm_write_unformatted,&   ! writes a full matrix to an open unit
00083             cp_fm_read_unformatted      ! reads a full matrix from an open unit
00084 
00085   PUBLIC :: cp_fm_indxg2p, &
00086             cp_fm_indxg2l, &
00087             cp_fm_indxl2g
00088 
00089   INTERFACE cp_fm_to_fm
00090     MODULE PROCEDURE cp_fm_to_fm_matrix, &  ! a full matrix
00091                      cp_fm_to_fm_columns    ! just a number of columns
00092   END INTERFACE
00093 
00094 ! *****************************************************************************
00109   TYPE cp_fm_type
00110 !    PRIVATE
00111      CHARACTER(LEN=60) :: name
00112      INTEGER :: id_nr, ref_count, print_count
00113      LOGICAL :: use_sp
00114      TYPE(cp_fm_struct_type), POINTER :: matrix_struct
00115      REAL(KIND = dp), DIMENSION(:,:), POINTER :: local_data
00116      REAL(KIND = sp), DIMENSION(:,:), POINTER :: local_data_sp
00117   END TYPE cp_fm_type
00118 
00119 ! *****************************************************************************
00126   TYPE cp_fm_p_type
00127      TYPE(cp_fm_type), POINTER :: matrix
00128   END TYPE cp_fm_p_type
00129 
00130 CONTAINS
00131 
00132 ! *****************************************************************************
00144   SUBROUTINE cp_fm_create(matrix,matrix_struct,name,use_sp,error)
00145     TYPE(cp_fm_type), POINTER                :: matrix
00146     TYPE(cp_fm_struct_type), POINTER         :: matrix_struct
00147     CHARACTER(len=*), INTENT(in), OPTIONAL   :: name
00148     LOGICAL, INTENT(in), OPTIONAL            :: use_sp
00149     TYPE(cp_error_type), INTENT(inout)       :: error
00150 
00151     CHARACTER(len=*), PARAMETER :: routineN = 'cp_fm_create', 
00152       routineP = moduleN//':'//routineN
00153 
00154     INTEGER                                  :: handle, ncol_local, npcol, 
00155                                                 nprow, nrow_local, stat
00156     LOGICAL                                  :: failure
00157     TYPE(cp_blacs_env_type), POINTER         :: context
00158 
00159     failure=.FALSE.
00160 
00161     CALL timeset(routineN,handle)
00162 
00163 #if defined(__parallel) && ! defined(__SCALAPACK)
00164      CALL cp_assert(.FALSE.,cp_failure_level,cp_assertion_failed,&
00165      routineP,"full matrixes need scalapack for parallel runs "//&
00166 CPSourceFileRef,&
00167      error)
00168 #endif
00169 
00170     CPPrecondition(ASSOCIATED(matrix_struct),cp_failure_level,routineP,error,failure)
00171     IF (.NOT.failure) THEN
00172        ALLOCATE(matrix,stat=stat)
00173        CPPostcondition(stat==0,cp_failure_level,routineP,error,failure)
00174     END IF
00175     IF (.NOT. failure) THEN
00176        context => matrix_struct%context
00177        matrix%matrix_struct => matrix_struct
00178        CALL cp_fm_struct_retain(matrix%matrix_struct,error=error)
00179        last_fm_id_nr=last_fm_id_nr+1
00180        matrix%id_nr=last_fm_id_nr
00181        matrix%ref_count=1
00182        matrix%print_count=0
00183 
00184        matrix%use_sp = .FALSE.
00185        IF(PRESENT(use_sp)) matrix%use_sp = use_sp
00186 
00187        nprow=context%num_pe(1)
00188        npcol=context%num_pe(2)
00189        NULLIFY(matrix%local_data)
00190        NULLIFY(matrix%local_data_sp)
00191 
00192        ! OK, we allocate here at least a 1 x 1 matrix
00193        ! this must (and is) compatible with the descinit call
00194        ! in cp_fm_struct
00195        nrow_local=matrix_struct%local_leading_dimension
00196        ncol_local=MAX(1,matrix_struct%ncol_locals(context%mepos(2)))
00197        IF(matrix%use_sp) THEN
00198           ALLOCATE(matrix%local_data_sp(nrow_local,ncol_local),stat=stat)
00199           CPPostcondition(stat==0,cp_failure_level,routineP,error,failure)
00200        ELSE
00201           ALLOCATE(matrix%local_data(nrow_local,ncol_local),stat=stat)
00202           CPPostcondition(stat==0,cp_failure_level,routineP,error,failure)
00203        ENDIF
00204 
00205        ! JVDV we should remove this, as it is up to the user to zero afterwards
00206        IF(matrix%use_sp) THEN
00207           CALL scopy(nrow_local*ncol_local,0.0_sp,0,matrix%local_data_sp,1)
00208        ELSE
00209           CALL dcopy(nrow_local*ncol_local,0.0_dp,0,matrix%local_data,1)
00210        ENDIF
00211 
00212        IF (PRESENT(name)) THEN
00213           matrix%name=name
00214        ELSE
00215           matrix%name='full matrix'//cp_to_string(matrix%id_nr)
00216        END IF
00217     END IF
00218 
00219     CALL timestop(handle)
00220 
00221   END SUBROUTINE cp_fm_create
00222 
00223 ! *****************************************************************************
00232   SUBROUTINE cp_fm_retain(matrix,error)
00233     TYPE(cp_fm_type), POINTER                :: matrix
00234     TYPE(cp_error_type), INTENT(inout)       :: error
00235 
00236     CHARACTER(len=*), PARAMETER :: routineN = 'cp_fm_retain', 
00237       routineP = moduleN//':'//routineN
00238 
00239     LOGICAL                                  :: failure
00240 
00241     failure=.FALSE.
00242 
00243     CPPrecondition(ASSOCIATED(matrix),cp_failure_level,routineP,error,failure)
00244     IF (.NOT. failure) THEN
00245        CPPrecondition(matrix%ref_count>0,cp_failure_level,routineP,error,failure)
00246        matrix%ref_count=matrix%ref_count+1
00247     END IF
00248 
00249   END SUBROUTINE cp_fm_retain
00250 
00251 ! *****************************************************************************
00260   SUBROUTINE cp_fm_release(matrix,error)
00261     TYPE(cp_fm_type), POINTER                :: matrix
00262     TYPE(cp_error_type), INTENT(inout)       :: error
00263 
00264     CHARACTER(len=*), PARAMETER :: routineN = 'cp_fm_release', 
00265       routineP = moduleN//':'//routineN
00266 
00267     INTEGER                                  :: handle, stat
00268     LOGICAL                                  :: failure
00269 
00270     failure=.FALSE.
00271 
00272     CALL timeset(routineN,handle)
00273 
00274     IF (ASSOCIATED(matrix)) THEN
00275        CPPrecondition(matrix%ref_count>0,cp_failure_level,routineP,error,failure)
00276        matrix%ref_count=matrix%ref_count-1
00277        IF (matrix%ref_count<1) THEN
00278           IF (ASSOCIATED(matrix%local_data)) THEN
00279              DEALLOCATE(matrix%local_data,stat=stat)
00280              CPPostcondition(stat==0,cp_warning_level,routineP,error,failure)
00281           END IF
00282           IF (ASSOCIATED(matrix%local_data_sp)) THEN
00283              DEALLOCATE(matrix%local_data_sp,stat=stat)
00284              CPPostcondition(stat==0,cp_warning_level,routineP,error,failure)
00285           END IF
00286           matrix%name=""
00287           CALL cp_fm_struct_release(matrix%matrix_struct,error=error)
00288           DEALLOCATE(matrix,stat=stat)
00289           CPPostcondition(stat==0,cp_warning_level,routineP,error,failure)
00290        END IF
00291     END IF
00292     NULLIFY(matrix)
00293 
00294     CALL timestop(handle)
00295 
00296   END SUBROUTINE cp_fm_release
00297 
00298 ! *****************************************************************************
00307   SUBROUTINE cp_fm_init_random(matrix,ncol,start_col,error)
00308     TYPE(cp_fm_type), POINTER                :: matrix
00309     INTEGER, INTENT(IN), OPTIONAL            :: ncol, start_col
00310     TYPE(cp_error_type), INTENT(inout)       :: error
00311 
00312     CHARACTER(len=*), PARAMETER :: routineN = 'cp_fm_init_random', 
00313       routineP = moduleN//':'//routineN
00314 
00315     INTEGER :: handle, icol_global, icol_local, irow_local, my_ncol, 
00316       my_start_col, ncol_global, ncol_local, nrow_global, nrow_local
00317     INTEGER, DIMENSION(:), POINTER           :: col_indices, row_indices
00318     LOGICAL                                  :: failure
00319     LOGICAL, SAVE                            :: FIRST = .TRUE.
00320     REAL(KIND=dp), ALLOCATABLE, DIMENSION(:) :: buff
00321     REAL(KIND=dp), DIMENSION(3, 2), SAVE     :: seed
00322     REAL(KIND=dp), DIMENSION(:, :), POINTER  :: local_data
00323     TYPE(rng_stream_type), POINTER           :: rng
00324 
00325     failure = .FALSE.
00326 
00327     CALL timeset(routineN,handle)
00328 
00329     ! first setup a parallel rng
00330     IF (FIRST) THEN
00331        seed(:,:)=RESHAPE((/1.0_dp,2.0_dp,3.0_dp,4.0_dp,5.0_dp,6.0_dp/),(/3,2/))
00332        FIRST=.FALSE.
00333     ENDIF
00334     ! guarantee same seed over all tasks
00335     CALL mp_bcast(seed,0,matrix%matrix_struct%para_env%group)
00336 
00337     NULLIFY(rng)
00338     CALL create_rng_stream(rng,"cp_fm_init_random_stream",distribution_type=UNIFORM, &
00339                            extended_precision=.TRUE.,seed=seed,error=error)
00340 
00341     CPPrecondition(.NOT.matrix%use_sp,cp_failure_level,routineP,error,failure)
00342 
00343     CALL cp_fm_get_info(matrix,nrow_global=nrow_global,ncol_global=ncol_global, &
00344                                nrow_local=nrow_local,ncol_local=ncol_local,&
00345                                local_data=local_data,&
00346                                row_indices=row_indices, col_indices=col_indices, error=error)
00347 
00348     my_start_col = 1
00349     IF (PRESENT(start_col)) my_start_col=start_col
00350     my_ncol = matrix%matrix_struct%ncol_global
00351     IF (PRESENT(ncol)) my_ncol=ncol
00352 
00353     CALL cp_assert(ncol_global>=(my_start_col+my_ncol-1),&
00354          cp_failure_level,cp_assertion_failed,routineP,&
00355          "ncol_global>=(my_start_col+my_ncol-1)",&
00356          error,failure)
00357 
00358     ALLOCATE(buff(nrow_global))
00359 
00360     ! each global row has its own substream, in order to reach the stream for the local col,
00361     ! we just reset to the next substream
00362     ! following this, we fill the full buff with random numbers, and pick those we need
00363     icol_global=0
00364     DO icol_local=1,ncol_local
00365        CPPrecondition(col_indices(icol_local)>icol_global,cp_failure_level,routineP,error,failure)
00366        DO
00367          CALL reset_to_next_rng_substream(rng,error=error)
00368          icol_global=icol_global+1
00369          IF (icol_global==col_indices(icol_local)) EXIT
00370        ENDDO
00371        CALL random_numbers(buff,rng,error)
00372        DO irow_local=1,nrow_local
00373           local_data(irow_local,icol_local)=buff(row_indices(irow_local))
00374        ENDDO
00375     ENDDO
00376 
00377     DEALLOCATE(buff)
00378 
00379     ! store seed before deletion (unclear if this is the proper seed)
00380     CALL get_rng_stream(rng,ig=seed,error=error)
00381     CALL delete_rng_stream(rng,error)
00382 
00383     CALL timestop(handle)
00384 
00385   END SUBROUTINE cp_fm_init_random
00386 
00387 ! *****************************************************************************
00394   SUBROUTINE cp_fm_set_all(matrix,alpha,beta,error)
00395 
00396     TYPE(cp_fm_type), POINTER                :: matrix
00397     REAL(KIND=dp), INTENT(IN)                :: alpha
00398     REAL(KIND=dp), INTENT(IN), OPTIONAL      :: beta
00399     TYPE(cp_error_type), INTENT(inout)       :: error
00400 
00401     CHARACTER(len=*), PARAMETER :: routineN = 'cp_fm_set_all', 
00402       routineP = moduleN//':'//routineN
00403 
00404     INTEGER                                  :: i, matrix_size, n
00405     LOGICAL                                  :: failure
00406 
00407 ! -------------------------------------------------------------------------
00408 
00409     failure = .FALSE.
00410 
00411     IF(matrix%use_sp) THEN
00412        matrix_size = SIZE(matrix%local_data_sp,1)*SIZE(matrix%local_data_sp,2)
00413        CALL scopy (matrix_size, REAL(alpha,sp), 0, matrix%local_data_sp, 1)
00414     ELSE
00415        matrix_size = SIZE(matrix%local_data,1)*SIZE(matrix%local_data,2)
00416        CALL dcopy (matrix_size, alpha, 0, matrix%local_data, 1)
00417     ENDIF
00418 
00419     IF (PRESENT(beta)) THEN
00420        CPPrecondition(.NOT.matrix%use_sp,cp_failure_level,routineP,error,failure)
00421       n = MIN(matrix%matrix_struct%nrow_global,matrix%matrix_struct%ncol_global)
00422       DO i=1,n
00423          CALL cp_fm_set_element(matrix,i,i,beta,error=error)
00424       END DO
00425     END IF
00426 
00427   END SUBROUTINE cp_fm_set_all
00428 
00429 ! *****************************************************************************
00432   SUBROUTINE cp_fm_get_diag(matrix,diag,error)
00433 
00434     IMPLICIT NONE
00435 
00436     ! arguments
00437     TYPE(cp_fm_type), POINTER                  :: matrix
00438     REAL(KIND = dp), DIMENSION(:), INTENT(OUT) :: diag
00439     TYPE(cp_error_type), INTENT(inout)         :: error
00440 
00441     CHARACTER(len=*), PARAMETER :: routineN = 'cp_fm_get_diag', 
00442       routineP = moduleN//':'//routineN
00443 
00444     ! locals
00445     INTEGER :: i,nrow_global
00446     LOGICAL :: failure
00447 
00448 #if defined(__SCALAPACK)
00449     INTEGER, DIMENSION(9) :: desca
00450     TYPE(cp_blacs_env_type), POINTER :: context
00451     INTEGER :: icol_local,ipcol,iprow,irow_local,mypcol,myprow,npcol,
00452                nprow
00453     REAL(KIND = dp), DIMENSION(:,:), POINTER :: a
00454     REAL(KIND = sp), DIMENSION(:,:), POINTER :: a_sp
00455 #endif
00456 
00457     failure = .FALSE.
00458 
00459     CALL cp_fm_get_info(matrix,nrow_global=nrow_global,error=error)
00460 
00461 #if defined(__SCALAPACK)
00462     diag=0.0_dp
00463     context => matrix%matrix_struct%context
00464     myprow=context%mepos(1)
00465     mypcol=context%mepos(2)
00466     nprow=context%num_pe(1)
00467     npcol=context%num_pe(2)
00468 
00469     a => matrix%local_data
00470     a_sp => matrix%local_data_sp
00471     desca(:) = matrix%matrix_struct%descriptor(:)
00472 
00473     DO i=1,nrow_global
00474        CALL infog2l(i,i,desca,nprow,npcol,myprow,mypcol,&
00475                     irow_local,icol_local,iprow,ipcol)
00476        IF ((iprow == myprow).AND.(ipcol == mypcol)) THEN
00477           IF(matrix%use_sp) THEN
00478              diag(i) = REAL(a_sp(irow_local,icol_local),dp)
00479           ELSE
00480              diag(i) = a(irow_local,icol_local)
00481           ENDIF
00482        ENDIF
00483     END DO
00484 #else
00485     DO i=1, nrow_global
00486        IF(matrix%use_sp) THEN
00487           diag(i)=REAL(matrix%local_data_sp(i,i),dp)
00488        ELSE
00489           diag(i)=matrix%local_data(i,i)
00490        ENDIF
00491     ENDDO
00492 #endif
00493     CALL mp_sum(diag,matrix%matrix_struct%para_env%group)
00494 
00495   END SUBROUTINE cp_fm_get_diag
00496 
00497 ! *****************************************************************************
00500   SUBROUTINE cp_fm_set_diag(matrix,diag,error)
00501 
00502     IMPLICIT NONE
00503 
00504     ! arguments
00505     TYPE(cp_fm_type), POINTER                  :: matrix
00506     REAL(KIND = dp), DIMENSION(:), INTENT(IN) :: diag
00507     TYPE(cp_error_type), INTENT(inout)         :: error
00508 
00509     CHARACTER(len=*), PARAMETER :: routineN = 'cp_fm_set_diag', 
00510       routineP = moduleN//':'//routineN
00511 
00512     ! locals
00513     INTEGER :: i,nrow_global
00514     LOGICAL :: failure
00515 #if defined(__SCALAPACK)
00516     INTEGER, DIMENSION(9) :: desca
00517     TYPE(cp_blacs_env_type), POINTER :: context
00518     INTEGER :: icol_local,ipcol,iprow,irow_local,mypcol,myprow,npcol,
00519                nprow
00520     REAL(KIND = dp), DIMENSION(:,:), POINTER :: a
00521 #endif
00522 
00523     failure = .FALSE.
00524 
00525     CALL cp_fm_get_info(matrix,nrow_global=nrow_global,error=error)
00526 
00527     CPPrecondition(.NOT.matrix%use_sp,cp_failure_level,routineP,error,failure)
00528 
00529 #if defined(__SCALAPACK)
00530     context => matrix%matrix_struct%context
00531     myprow=context%mepos(1)
00532     mypcol=context%mepos(2)
00533     nprow=context%num_pe(1)
00534     npcol=context%num_pe(2)
00535 
00536     a => matrix%local_data
00537     desca(:) = matrix%matrix_struct%descriptor(:)
00538 
00539     DO i=1,nrow_global
00540        CALL infog2l(i,i,desca,nprow,npcol,myprow,mypcol,&
00541                     irow_local,icol_local,iprow,ipcol)
00542        IF ((iprow == myprow).AND.(ipcol == mypcol)) THEN
00543          a(irow_local,icol_local) = diag(i)
00544        ENDIF
00545     END DO
00546 #else
00547     DO i=1, nrow_global
00548        matrix%local_data(i,i)=diag(i)
00549     ENDDO
00550 #endif
00551 
00552   END SUBROUTINE cp_fm_set_diag
00553 
00554 
00555 
00556 
00557 ! *****************************************************************************
00574   SUBROUTINE cp_fm_get_element(matrix,irow_global,icol_global,alpha,local)
00575 
00576     IMPLICIT NONE
00577 
00578     ! arguments
00579     TYPE(cp_fm_type), POINTER          :: matrix
00580     REAL(KIND = dp), INTENT(OUT)                     :: alpha
00581     INTEGER, INTENT(IN)                       :: icol_global,
00582                                                  irow_global
00583     LOGICAL, INTENT(OUT), OPTIONAL            :: local
00584     LOGICAL :: failure
00585 
00586     ! locals
00587 #if defined(__SCALAPACK)
00588     INTEGER, DIMENSION(9) :: desca
00589     TYPE(cp_blacs_env_type), POINTER :: context
00590     INTEGER :: icol_local,ipcol,iprow,irow_local,mypcol,myprow,npcol,
00591                nprow
00592     REAL(KIND = dp), DIMENSION(:,:), POINTER :: a
00593 #endif
00594 
00595     failure = .FALSE.
00596 
00597 #if defined(__SCALAPACK)
00598     context => matrix%matrix_struct%context
00599     myprow=context%mepos(1)
00600     mypcol=context%mepos(2)
00601     nprow=context%num_pe(1)
00602     npcol=context%num_pe(2)
00603 
00604     a => matrix%local_data
00605     desca(:) = matrix%matrix_struct%descriptor(:)
00606 
00607     CALL infog2l(irow_global,icol_global,desca,nprow,npcol,myprow,mypcol,&
00608                  irow_local,icol_local,iprow,ipcol)
00609 
00610     IF ((iprow == myprow).AND.(ipcol == mypcol)) THEN
00611       alpha = a(irow_local,icol_local)
00612       CALL cp_blacs_dgebs2d(context%group, 'All', ' ', 1, 1, alpha , 1 )
00613       IF (PRESENT(local)) local=.TRUE.
00614     ELSE
00615       CALL cp_blacs_dgebr2d(context%group, 'All', ' ', 1, 1, alpha , 1 , iprow, ipcol )
00616       IF (PRESENT(local)) local=.FALSE.
00617     END IF
00618 
00619 #else
00620     IF (PRESENT(local)) local=.TRUE.
00621     alpha = matrix%local_data(irow_global,icol_global)
00622 #endif
00623 
00624   END SUBROUTINE cp_fm_get_element
00625 
00626 ! *****************************************************************************
00632   SUBROUTINE cp_fm_set_element(matrix,irow_global,icol_global,alpha,error)
00633     TYPE(cp_fm_type), POINTER                :: matrix
00634     INTEGER, INTENT(IN)                      :: irow_global, icol_global
00635     REAL(KIND=dp), INTENT(IN)                :: alpha
00636     TYPE(cp_error_type), INTENT(inout)       :: error
00637 
00638     CHARACTER(len=*), PARAMETER :: routineN = 'cp_fm_set_element', 
00639       routineP = moduleN//':'//routineN
00640 
00641     INTEGER                                  :: icol_local, ipcol, iprow, 
00642                                                 irow_local, mypcol, myprow, 
00643                                                 npcol, nprow
00644     INTEGER, DIMENSION(9)                    :: desca
00645     LOGICAL                                  :: failure
00646     REAL(KIND=dp), DIMENSION(:, :), POINTER  :: a
00647     TYPE(cp_blacs_env_type), POINTER         :: context
00648 
00649     failure = .FALSE.
00650 
00651     context => matrix%matrix_struct%context
00652     myprow=context%mepos(1)
00653     mypcol=context%mepos(2)
00654     nprow=context%num_pe(1)
00655     npcol=context%num_pe(2)
00656 
00657     CPPrecondition(.NOT.matrix%use_sp,cp_failure_level,routineP,error,failure)
00658 
00659 #if defined(__SCALAPACK)
00660 
00661     a => matrix%local_data
00662 
00663     desca(:) = matrix%matrix_struct%descriptor(:)
00664 
00665     CALL infog2l(irow_global,icol_global,desca,nprow,npcol,myprow,mypcol,&
00666          irow_local,icol_local,iprow,ipcol)
00667 
00668     IF ((iprow == myprow).AND.(ipcol == mypcol)) THEN
00669        a(irow_local,icol_local) = alpha
00670     END IF
00671 
00672 #else
00673 
00674     matrix%local_data(irow_global,icol_global) = alpha
00675 
00676 #endif
00677   END SUBROUTINE cp_fm_set_element
00678 
00679 ! *****************************************************************************
00706   SUBROUTINE cp_fm_set_submatrix(fm,new_values,start_row,&
00707        start_col, n_rows, n_cols, alpha, beta, transpose, error)
00708     TYPE(cp_fm_type), POINTER                :: fm
00709     REAL(KIND=dp), DIMENSION(:, :), 
00710       INTENT(in)                             :: new_values
00711     INTEGER, INTENT(in), OPTIONAL            :: start_row, start_col, n_rows, 
00712                                                 n_cols
00713     REAL(KIND=dp), INTENT(in), OPTIONAL      :: alpha, beta
00714     LOGICAL, INTENT(in), OPTIONAL            :: transpose
00715     TYPE(cp_error_type), INTENT(inout)       :: error
00716 
00717     CHARACTER(len=*), PARAMETER :: routineN = 'cp_fm_set_submatrix', 
00718       routineP = moduleN//':'//routineN
00719 
00720     INTEGER :: i, i0, j, j0, ncol, ncol_block, ncol_global, ncol_local, nrow, 
00721       nrow_block, nrow_global, nrow_local, this_col, this_row
00722     INTEGER, DIMENSION(:), POINTER           :: col_indices, row_indices
00723     LOGICAL                                  :: failure, tr_a
00724     REAL(KIND=dp)                            :: al, be
00725     REAL(KIND=dp), DIMENSION(:, :), POINTER  :: full_block
00726 
00727     al=1.0_dp; be=0.0_dp; i0=1; j0=1; tr_a=.FALSE.
00728     ! can be called too many times, making it a bit useless
00729     ! CALL timeset(routineN//','//moduleN,handle)
00730 
00731     failure=.FALSE.
00732 
00733     CPPrecondition(ASSOCIATED(fm),cp_failure_level,routineP,error,failure)
00734     IF (.not.failure) THEN
00735        CPPrecondition(fm%ref_count>0,cp_failure_level,routineP,error,failure)
00736     END IF
00737 
00738     CPPrecondition(.NOT.fm%use_sp,cp_failure_level,routineP,error,failure)
00739 
00740     IF (.not.failure) THEN
00741        IF (PRESENT(alpha)) al=alpha
00742        IF (PRESENT(beta)) be=beta
00743        IF (PRESENT(start_row)) i0=start_row
00744        IF (PRESENT(start_col)) j0=start_col
00745        IF (PRESENT(transpose)) tr_a=transpose
00746        IF (tr_a) THEN
00747           nrow=SIZE(new_values,2)
00748           ncol=SIZE(new_values,1)
00749        ELSE
00750           nrow=SIZE(new_values,1)
00751           ncol=SIZE(new_values,2)
00752        END IF
00753        IF (PRESENT(n_rows)) nrow=n_rows
00754        IF (PRESENT(n_cols)) ncol=n_cols
00755 
00756        full_block => fm%local_data
00757 
00758        CALL cp_fm_get_info(matrix=fm,&
00759             nrow_global=nrow_global,ncol_global=ncol_global,&
00760             nrow_block =nrow_block ,ncol_block =ncol_block ,&
00761             nrow_local =nrow_local ,ncol_local =ncol_local ,&
00762             row_indices=row_indices,col_indices=col_indices,error=error)
00763 
00764        IF (al==1.0.AND.be==0.0) THEN
00765           DO j=1,ncol_local
00766              this_col=col_indices(j)-j0+1
00767              IF (this_col.GE.1 .AND. this_col.LE.ncol) THEN
00768                 IF (tr_a) THEN
00769                    IF (i0==1.AND.nrow_global==nrow) THEN
00770                       DO i=1,nrow_local
00771                          full_block(i,j)=new_values(this_col,row_indices(i))
00772                       END DO
00773                    ELSE
00774                       DO i=1,nrow_local
00775                          this_row=row_indices(i)-i0+1
00776                          IF (this_row>=1 .AND. this_row<=nrow) THEN
00777                             full_block(i,j)=new_values(this_col,this_row)
00778                          END IF
00779                       END DO
00780                    END IF
00781                 ELSE
00782                    IF (i0==1.AND.nrow_global==nrow) THEN
00783                       DO i=1,nrow_local
00784                          full_block(i,j)=new_values(row_indices(i),this_col)
00785                       END DO
00786                    ELSE
00787                       DO i=1,nrow_local
00788                          this_row=row_indices(i)-i0+1
00789                          IF (this_row>=1 .AND. this_row<=nrow) THEN
00790                             full_block(i,j)=new_values(this_row,this_col)
00791                          END IF
00792                       END DO
00793                    END IF
00794                 END IF
00795              END IF
00796           END DO
00797        ELSE
00798           DO j=1,ncol_local
00799              this_col=col_indices(j)-j0+1
00800              IF (this_col.GE.1 .AND. this_col.LE.ncol) THEN
00801                 IF (tr_a) THEN
00802                    DO i=1,nrow_local
00803                       this_row=row_indices(i)-i0+1
00804                       IF (this_row>=1 .AND. this_row<=nrow) THEN
00805                          full_block(i,j)=al*new_values(this_col,this_row)+&
00806                               be*full_block(i,j)
00807                       END IF
00808                    END DO
00809                 ELSE
00810                    DO i=1,nrow_local
00811                       this_row=row_indices(i)-i0+1
00812                       IF (this_row>=1 .AND. this_row<=nrow) THEN
00813                          full_block(i,j)=al*new_values(this_row,this_col)+&
00814                               be*full_block(i,j)
00815                       END IF
00816                    END DO
00817                 END IF
00818              END IF
00819           END DO
00820        END IF
00821     END IF
00822 
00823     ! CALL timestop(handle)
00824 
00825   END SUBROUTINE cp_fm_set_submatrix
00826 
00827 ! *****************************************************************************
00854   SUBROUTINE cp_fm_get_submatrix(fm,target_m, start_row,&
00855        start_col, n_rows, n_cols, transpose, error)
00856     TYPE(cp_fm_type), POINTER                :: fm
00857     REAL(KIND=dp), DIMENSION(:, :), 
00858       INTENT(out)                            :: target_m
00859     INTEGER, INTENT(in), OPTIONAL            :: start_row, start_col, n_rows, 
00860                                                 n_cols
00861     LOGICAL, INTENT(in), OPTIONAL            :: transpose
00862     TYPE(cp_error_type), INTENT(inout)       :: error
00863 
00864     CHARACTER(len=*), PARAMETER :: routineN = 'cp_fm_get_submatrix', 
00865       routineP = moduleN//':'//routineN
00866 
00867     INTEGER :: handle, i, i0, j, j0, ncol, ncol_global, ncol_local, nrow, 
00868       nrow_global, nrow_local, this_col, this_row
00869     INTEGER, DIMENSION(:), POINTER           :: col_indices, row_indices
00870     LOGICAL                                  :: failure, tr_a
00871     REAL(KIND=dp), DIMENSION(:, :), POINTER  :: full_block
00872     TYPE(cp_para_env_type), POINTER          :: para_env
00873 
00874     CALL timeset(routineN,handle)
00875 
00876     i0=1; j0=1; tr_a=.FALSE.
00877     failure=.FALSE.
00878 
00879     CPPrecondition(ASSOCIATED(fm),cp_failure_level,routineP,error,failure)
00880     IF (.not.failure) THEN
00881        CPPrecondition(fm%ref_count>0,cp_failure_level,routineP,error,failure)
00882     END IF
00883 
00884     CPPrecondition(.NOT.fm%use_sp,cp_failure_level,routineP,error,failure)
00885 
00886     IF (.not.failure) THEN
00887        IF (PRESENT(start_row)) i0=start_row
00888        IF (PRESENT(start_col)) j0=start_col
00889        IF (PRESENT(transpose)) tr_a=transpose
00890        IF (tr_a) THEN
00891           nrow=SIZE(target_m,2)
00892           ncol=SIZE(target_m,1)
00893        ELSE
00894           nrow=SIZE(target_m,1)
00895           ncol=SIZE(target_m,2)
00896        END IF
00897        IF (PRESENT(n_rows)) nrow=n_rows
00898        IF (PRESENT(n_cols)) ncol=n_cols
00899 
00900        para_env => fm%matrix_struct%para_env
00901 
00902        full_block => fm%local_data
00903 #if defined(__SCALAPACK)
00904        ! zero-out whole target_m
00905        IF (SIZE(target_m,1)*SIZE(target_m,2)/=0) THEN
00906           CALL dcopy(SIZE(target_m,1)*SIZE(target_m,2),0.0_dp,0,target_m(1,1),1)
00907        END IF
00908 #endif
00909 
00910        CALL cp_fm_get_info(matrix=fm, &
00911             nrow_global=nrow_global,ncol_global=ncol_global, &
00912             nrow_local =nrow_local ,ncol_local =ncol_local , &
00913             row_indices=row_indices,col_indices=col_indices,error=error)
00914 
00915        DO j=1,ncol_local
00916           this_col=col_indices(j)-j0+1
00917           IF (this_col.GE.1 .AND. this_col.LE.ncol) THEN
00918              IF (tr_a) THEN
00919                 IF (i0==1.AND.nrow_global==nrow) THEN
00920                    DO i=1,nrow_local
00921                       target_m(this_col,row_indices(i))=full_block(i,j)
00922                    END DO
00923                 ELSE
00924                    DO i=1,nrow_local
00925                       this_row=row_indices(i)-i0+1
00926                       IF (this_row>=1 .AND. this_row<=nrow) THEN
00927                          target_m(this_col,this_row)=full_block(i,j)
00928                       END IF
00929                    END DO
00930                 END IF
00931              ELSE
00932                 IF (i0==1.AND.nrow_global==nrow) THEN
00933                    DO i=1,nrow_local
00934                       target_m(row_indices(i),this_col)=full_block(i,j)
00935                    END DO
00936                 ELSE
00937                    DO i=1,nrow_local
00938                       this_row=row_indices(i)-i0+1
00939                       IF (this_row>=1 .AND. this_row<=nrow) THEN
00940                          target_m(this_row,this_col)=full_block(i,j)
00941                       END IF
00942                    END DO
00943                 END IF
00944              END IF
00945           END IF
00946        END DO
00947 
00948        CALL mp_sum(target_m,para_env%group)
00949     END IF
00950 
00951     CALL timestop(handle)
00952 
00953   END SUBROUTINE cp_fm_get_submatrix
00954 
00955 ! *****************************************************************************
00962   SUBROUTINE cp_fm_get_info(matrix,name,nrow_global,ncol_global,&
00963        nrow_block,ncol_block,nrow_local,ncol_local,&
00964        row_indices,col_indices,local_data,context,&
00965        nrow_locals, ncol_locals, matrix_struct,para_env,error)
00966 
00967     TYPE(cp_fm_type), POINTER                  :: matrix
00968     CHARACTER(LEN=*), OPTIONAL, INTENT(OUT)    :: name
00969     INTEGER, OPTIONAL, INTENT(OUT)             :: ncol_block,ncol_global,
00970                                                   nrow_block,nrow_global,
00971                                                   nrow_local,ncol_local
00972     INTEGER, OPTIONAL, DIMENSION(:), POINTER   :: row_indices,col_indices,
00973                                                   nrow_locals, ncol_locals
00974     TYPE(cp_para_env_type), POINTER, OPTIONAL  :: para_env
00975     TYPE(cp_blacs_env_type), POINTER, OPTIONAL :: context
00976     TYPE(cp_fm_struct_type),POINTER,OPTIONAL   :: matrix_struct
00977     TYPE(cp_error_type),INTENT(inout):: error
00978     REAL(KIND = dp), DIMENSION(:,:),POINTER, OPTIONAL :: local_data
00979 
00980     CHARACTER(len=*), PARAMETER :: routineN='cp_fm_get_info',
00981          routineP=moduleN//':'//routineN
00982     INTEGER i,nprow,npcol,myprow,mypcol, stat
00983     TYPE(cp_blacs_env_type), POINTER :: ctxt
00984     LOGICAL :: failure
00985 #if defined(__SCALAPACK)
00986     INTEGER , EXTERNAL :: indxl2g
00987 #endif
00988 
00989     IF (PRESENT(name)) name = matrix%name
00990     IF (PRESENT(context)) context => matrix%matrix_struct%context
00991     IF (PRESENT(matrix_struct)) matrix_struct => matrix%matrix_struct
00992     IF (PRESENT(local_data)) local_data => matrix%local_data ! not hiding things anymore :-(
00993     IF (PRESENT(row_indices)) THEN
00994        CALL cp_fm_struct_get(matrix%matrix_struct, row_indices=row_indices,&
00995             error=error)
00996     ENDIF
00997     IF (PRESENT(col_indices)) THEN
00998        CALL cp_fm_struct_get(matrix%matrix_struct, col_indices=col_indices,&
00999             error=error)
01000     ENDIF
01001     IF (PRESENT(nrow_locals)) THEN
01002        nrow_locals => matrix%matrix_struct%nrow_locals
01003     END IF
01004     IF (PRESENT(ncol_locals)) THEN
01005        nrow_locals => matrix%matrix_struct%ncol_locals
01006     END IF
01007 
01008     CALL cp_fm_struct_get(matrix%matrix_struct, nrow_local=nrow_local,&
01009           ncol_local=ncol_local, nrow_global=nrow_global,&
01010           ncol_global=ncol_global, nrow_block=nrow_block, &
01011           ncol_block=ncol_block, error=error)
01012 
01013     IF (PRESENT(para_env)) para_env => matrix%matrix_struct%para_env
01014 
01015   END SUBROUTINE cp_fm_get_info
01016 
01017 ! *****************************************************************************
01021   SUBROUTINE cp_fm_maxabsval(matrix,a_max,ir_max,ic_max,error)
01022     TYPE(cp_fm_type), POINTER                :: matrix
01023     REAL(KIND=dp), INTENT(OUT)               :: a_max
01024     INTEGER, INTENT(OUT), OPTIONAL           :: ir_max, ic_max
01025     TYPE(cp_error_type), INTENT(inout)       :: error
01026 
01027     CHARACTER(len=*), PARAMETER :: routineN = 'cp_fm_maxabsval', 
01028       routineP = moduleN//':'//routineN
01029 
01030     INTEGER                                  :: handle, i, ic_max_local, 
01031                                                 ir_max_local, istat, j, 
01032                                                 mepos, ncol_local, 
01033                                                 nrow_local, num_pe
01034     INTEGER, ALLOCATABLE, DIMENSION(:)       :: ic_max_vec, ir_max_vec
01035     INTEGER, DIMENSION(:), POINTER           :: col_indices, row_indices
01036     LOGICAL                                  :: failure
01037     REAL(dp)                                 :: my_max
01038     REAL(dp), ALLOCATABLE, DIMENSION(:)      :: a_max_vec
01039     REAL(KIND=dp), DIMENSION(:, :), POINTER  :: my_block
01040     REAL(KIND=sp), DIMENSION(:, :), POINTER  :: my_block_sp
01041 
01042     CALL timeset(routineN,handle)
01043 
01044     failure = .FALSE.
01045     my_block => matrix%local_data
01046     my_block_sp => matrix%local_data_sp
01047 
01048     CALL cp_fm_get_info(matrix, nrow_local=nrow_local, ncol_local=ncol_local,&
01049                         row_indices=row_indices,col_indices=col_indices,error=error)
01050 
01051     IF(matrix%use_sp) THEN
01052        a_max = REAL(MAXVAL(ABS(my_block_sp(1:nrow_local,1:ncol_local))),dp)
01053     ELSE
01054        a_max = MAXVAL(ABS(my_block(1:nrow_local,1:ncol_local)))
01055     ENDIF
01056 
01057     IF(PRESENT(ir_max)) THEN
01058       num_pe = matrix%matrix_struct%para_env%num_pe
01059       mepos = matrix%matrix_struct%para_env%mepos
01060       ALLOCATE(ir_max_vec(0:num_pe-1),STAT=istat)
01061       CPPostcondition(istat==0,cp_failure_level,routineP,error,failure)
01062       ir_max_vec(0:num_pe-1) = 0
01063       ALLOCATE(ic_max_vec(0:num_pe-1),STAT=istat)
01064       CPPostcondition(istat==0,cp_failure_level,routineP,error,failure)
01065       ic_max_vec(0:num_pe-1) = 0
01066       ALLOCATE(a_max_vec(0:num_pe-1),STAT=istat)
01067       CPPostcondition(istat==0,cp_failure_level,routineP,error,failure)
01068       a_max_vec(0:num_pe-1) = 0.0_dp
01069       my_max = 0.0_dp
01070 
01071       IF((ncol_local > 0) .AND. (nrow_local > 0)) THEN
01072         DO i = 1,ncol_local
01073           DO j = 1,nrow_local
01074              IF(matrix%use_sp) THEN
01075                 IF( ABS(my_block_sp(j,i)) > my_max) THEN
01076                    my_max = REAL(my_block_sp(j,i),dp)
01077                    ir_max_local = j
01078                    ic_max_local = i
01079                 ENDIF
01080              ELSE
01081                 IF( ABS(my_block(j,i)) > my_max) THEN
01082                    my_max =  my_block(j,i)
01083                    ir_max_local = j
01084                    ic_max_local = i
01085                 ENDIF
01086             END IF
01087           END DO
01088         END DO
01089 
01090         a_max_vec(mepos) = my_max
01091         ir_max_vec(mepos) = row_indices(ir_max_local)
01092         ic_max_vec(mepos) = col_indices(ic_max_local)
01093 
01094       END IF
01095 
01096       CALL mp_sum(a_max_vec,matrix%matrix_struct%para_env%group)
01097       CALL mp_sum(ir_max_vec,matrix%matrix_struct%para_env%group)
01098       CALL mp_sum(ic_max_vec,matrix%matrix_struct%para_env%group)
01099 
01100       my_max = 0.0_dp
01101       DO i = 0,num_pe-1
01102         IF(a_max_vec(i)>my_max) THEN
01103            ir_max = ir_max_vec(i)
01104            ic_max = ic_max_vec(i)
01105         END IF
01106       END DO
01107 
01108       DEALLOCATE(ir_max_vec,ic_max_vec,a_max_vec, STAT = istat)
01109       CPPostcondition(istat==0,cp_failure_level,routineP,error,failure)
01110       CPPostcondition(ic_max>0,cp_failure_level,routineP,error,failure)
01111       CPPostcondition(ir_max>0,cp_failure_level,routineP,error,failure)
01112 
01113     END IF
01114 
01115     CALL mp_max(a_max,matrix%matrix_struct%para_env%group)
01116 
01117     CALL timestop(handle)
01118 
01119   END SUBROUTINE cp_fm_maxabsval
01120 
01121 ! *****************************************************************************
01130   SUBROUTINE cp_fm_maxabsrownorm(matrix, a_max,error)
01131     TYPE(cp_fm_type), POINTER                :: matrix
01132     REAL(KIND=dp), INTENT(OUT)               :: a_max
01133     TYPE(cp_error_type), INTENT(inout)       :: error
01134 
01135     CHARACTER(len=*), PARAMETER :: routineN = 'cp_fm_maxabsrownorm', 
01136       routineP = moduleN//':'//routineN
01137 
01138     INTEGER                                  :: handle, i, j, ncol_local, 
01139                                                 nrow_global, nrow_local
01140     INTEGER, DIMENSION(:), POINTER           :: row_indices
01141     LOGICAL                                  :: failure
01142     REAL(KIND=dp), ALLOCATABLE, DIMENSION(:) :: values
01143     REAL(KIND=dp), DIMENSION(:, :), POINTER  :: my_block
01144 
01145     CALL timeset(routineN,handle)
01146 
01147     my_block => matrix%local_data
01148 
01149     CPPrecondition(.NOT.matrix%use_sp,cp_failure_level,routineP,error,failure)
01150 
01151     CALL cp_fm_get_info(matrix, row_indices=row_indices, nrow_global=nrow_global, &
01152                                 nrow_local=nrow_local, ncol_local=ncol_local,error=error)
01153 
01154     ! the efficiency could be improved by making use of the row-col distribution of scalapack
01155     ALLOCATE(values(nrow_global))
01156     values=0.0_dp
01157     DO j=1,ncol_local
01158        DO i=1,nrow_local
01159           values(row_indices(i))=values(row_indices(i))+ABS(my_block(i,j))
01160        ENDDO
01161     ENDDO
01162     CALL mp_sum(values,matrix%matrix_struct%para_env%group)
01163     a_max = MAXVAL(values)
01164     DEALLOCATE(values)
01165 
01166     CALL timestop(handle)
01167   END SUBROUTINE
01168 
01169 ! *****************************************************************************
01172   SUBROUTINE cp_fm_vectorsnorm(matrix, norm_array, error)
01173     TYPE(cp_fm_type), POINTER                :: matrix
01174     REAL(KIND=dp), DIMENSION(:), INTENT(OUT) :: norm_array
01175     TYPE(cp_error_type), INTENT(inout)       :: error
01176 
01177     CHARACTER(len=*), PARAMETER :: routineN = 'cp_fm_vectorsnorm', 
01178       routineP = moduleN//':'//routineN
01179 
01180     INTEGER                                  :: handle, i, j, ncol_global, 
01181                                                 ncol_local, nrow_local
01182     INTEGER, DIMENSION(:), POINTER           :: col_indices
01183     LOGICAL                                  :: failure
01184     REAL(KIND=dp), DIMENSION(:, :), POINTER  :: my_block
01185 
01186     failure = .FALSE.
01187 
01188     CALL timeset(routineN,handle)
01189 
01190     my_block => matrix%local_data
01191 
01192     CPPrecondition(.NOT.matrix%use_sp,cp_failure_level,routineP,error,failure)
01193 
01194     CALL cp_fm_get_info(matrix, col_indices=col_indices, ncol_global=ncol_global, &
01195                                 nrow_local=nrow_local, ncol_local=ncol_local, error=error)
01196 
01197     ! the efficiency could be improved by making use of the row-col distribution of scalapack
01198     norm_array=0.0_dp
01199     DO j=1,ncol_local
01200        DO i=1,nrow_local
01201           norm_array(col_indices(j))=norm_array(col_indices(j))+my_block(i,j)*my_block(i,j)
01202        ENDDO
01203     ENDDO
01204     CALL mp_sum(norm_array,matrix%matrix_struct%para_env%group)
01205     norm_array = SQRT(norm_array)
01206 
01207     CALL timestop(handle)
01208   END SUBROUTINE
01209 
01210 ! *****************************************************************************
01229    SUBROUTINE cp_fm_write(matrix, unit_nr, long_description, local, error)
01230     TYPE(cp_fm_type), POINTER                :: matrix
01231     INTEGER, INTENT(in)                      :: unit_nr
01232     LOGICAL, INTENT(in), OPTIONAL            :: long_description, local
01233     TYPE(cp_error_type), INTENT(inout)       :: error
01234 
01235     CHARACTER(len=*), PARAMETER :: routineN = 'cp_fm_write', 
01236       routineP = moduleN//':'//routineN
01237 
01238     CHARACTER(len=default_path_length)       :: filename
01239     INTEGER                                  :: data_unit, desc_unit, iostat
01240     LOGICAL                                  :: exists, failure, loc, 
01241                                                 my_long_description
01242     TYPE(cp_logger_type), POINTER            :: logger
01243     TYPE(cp_para_env_type), POINTER          :: para_env
01244 
01245      failure=.FALSE.; loc=.FALSE.; my_long_description=.FALSE.
01246      filename=' '
01247 
01248      IF (PRESENT(local)) loc=local
01249      IF (PRESENT(long_description)) my_long_description=long_description
01250      logger=>cp_error_get_logger(error)
01251      para_env=>logger%para_env
01252      IF (ASSOCIATED(matrix)) THEN
01253         data_unit=unit_nr
01254         desc_unit=unit_nr
01255         IF (.not.loc.and.long_description) THEN
01256            matrix%print_count=matrix%print_count+1
01257            IF (para_env%mepos==para_env%source) THEN
01258               filename=TRIM(matrix%name)//'-'//&
01259                    ADJUSTL(cp_to_string(matrix%print_count))
01260               CALL compress(filename,full=.TRUE.)
01261               WRITE (unit=desc_unit,&
01262                    fmt="(' <cp_fm_types>:{ id_nr=',i10,' written to files ',a,'*')",&
01263                    iostat=iostat) matrix%id_nr, TRIM(filename)
01264            END IF
01265            CALL cp_logger_generate_filename(logger=logger,res=filename,&
01266                 root=TRIM(matrix%name)//"-"//&
01267                 TRIM(ADJUSTL(cp_to_string(matrix%print_count))),&
01268                 postfix=".desc",local=.TRUE.)
01269            INQUIRE (FILE=TRIM(filename),EXIST=exists)
01270            CALL cp_assert(.NOT.exists,cp_warning_level,cp_assertion_failed,&
01271                 routineP,"file "//TRIM(filename)//" exists, overwriting",error)
01272            CALL open_file(TRIM(filename),file_status="unknown",&
01273                 file_action="WRITE",&
01274                 unit_number=desc_unit)
01275 
01276            CALL cp_logger_generate_filename(logger=logger,res=filename,&
01277              root=TRIM(matrix%name)//"-"//&
01278              TRIM(ADJUSTL(cp_to_string(matrix%print_count))),&
01279              postfix=".dat",local=.TRUE.)
01280            INQUIRE (FILE=TRIM(filename),EXIST=exists)
01281            CALL cp_assert(.NOT.exists,cp_warning_level,cp_assertion_failed,&
01282                 routineP,"file "//TRIM(filename)//" exists, overwriting",error)
01283            CALL open_file(TRIM(filename),file_status="unknown",&
01284                 file_action="WRITE",&
01285                 unit_number=data_unit)
01286         END IF
01287 
01288         IF (loc .OR. para_env%mepos==para_env%source.or.long_description) THEN
01289            WRITE (unit=desc_unit,&
01290                 fmt="(' <cp_fm_types>:{ id_nr=',i10,' ref_count=',i10,',')",&
01291                 iostat=iostat) matrix%id_nr, matrix%ref_count
01292            WRITE (unit=desc_unit,&
01293                 fmt="(' name=',a,',')",&
01294                 iostat=iostat) matrix%name
01295            WRITE(unit=desc_unit,fmt="(a)",iostat=iostat) " matrix_structure="
01296            CALL cp_fm_struct_write(matrix%matrix_struct,unit_nr=desc_unit,&
01297                 long_description=my_long_description,error=error)
01298            WRITE(unit=desc_unit,fmt="(a)",iostat=iostat) " ,"
01299            IF (my_long_description) THEN
01300               WRITE(unit=desc_unit,fmt="(a)",iostat=iostat) " local_data=("
01301               CALL cp_2d_r_write(matrix%local_data,unit_nr=data_unit,error=error)
01302               IF (desc_unit/=data_unit) THEN
01303                  WRITE(unit=desc_unit,fmt="(a)",iostat=iostat) &
01304                       "  written to file "//TRIM(filename)
01305               END IF
01306            ELSE IF (loc) THEN
01307               IF (ASSOCIATED(matrix%local_data)) THEN
01308                  WRITE (unit=desc_unit,&
01309                       fmt="(' local_data=(REAL(wp, DIMENSION(',i8,', ',i8,'))')",&
01310                       iostat=iostat) SIZE(matrix%local_data,1), SIZE(matrix%local_data,2)
01311               END IF
01312            END IF
01313            WRITE (unit=desc_unit,fmt="(a)",iostat=iostat) " }"
01314         END IF
01315 
01316         IF (.not.loc.and.long_description) THEN
01317            CALL close_file(desc_unit)
01318            CALL close_file(data_unit)
01319         END IF
01320      ELSE
01321         IF (loc .OR. para_env%mepos==para_env%source) THEN
01322            WRITE (unit=unit_nr,fmt="(a)",iostat=iostat) " <cp_fm_types>:*null*"
01323         END IF
01324      END IF
01325    END SUBROUTINE cp_fm_write
01326 
01327 ! *****************************************************************************
01332   SUBROUTINE cp_fm_to_fm_matrix(source,destination,error)
01333 
01334     TYPE(cp_fm_type), POINTER                :: source, destination
01335     TYPE(cp_error_type), INTENT(inout)       :: error
01336 
01337     CHARACTER(LEN=*), PARAMETER :: routineN = 'cp_fm_to_fm_matrix', 
01338       routineP = moduleN//':'//routineN
01339 
01340     INTEGER                                  :: handle, npcol, nprow
01341     LOGICAL                                  :: failure
01342 
01343     CALL timeset(routineN,handle)
01344 
01345     failure=.TRUE.
01346 
01347     nprow = source%matrix_struct%context%num_pe(1)
01348     npcol = source%matrix_struct%context%num_pe(2)
01349 
01350     IF ((.NOT. cp2k_is_parallel).OR.&
01351          cp_fm_struct_equivalent(source%matrix_struct,&
01352          destination%matrix_struct,error=error)) THEN
01353        IF(source%use_sp.AND.destination%use_sp) THEN
01354           CALL cp_assert(SIZE(source%local_data_sp,1)==SIZE(destination%local_data_sp,1).AND.&
01355                SIZE(source%local_data_sp,2)==SIZE(destination%local_data_sp,2),cp_failure_level,&
01356                cp_assertion_failed, routineP,&
01357                "Cannot copy full matrix <"//TRIM(source%name)//&
01358                "> to full matrix <"//TRIM(destination%name)//&
01359                ">. The local_data blocks have different sizes."//&
01360 CPSourceFileRef,&
01361                error=error)
01362           CALL scopy(SIZE(source%local_data_sp,1)*SIZE(source%local_data_sp,2),&
01363                source%local_data_sp(1,1),1,destination%local_data_sp(1,1),1)
01364        ELSEIF(source%use_sp.AND..NOT.destination%use_sp) THEN
01365           CALL cp_assert(SIZE(source%local_data_sp,1)==SIZE(destination%local_data,1).AND.&
01366                SIZE(source%local_data_sp,2)==SIZE(destination%local_data,2),cp_failure_level,&
01367                cp_assertion_failed, routineP,&
01368                "Cannot copy full matrix <"//TRIM(source%name)//&
01369                "> to full matrix <"//TRIM(destination%name)//&
01370                ">. The local_data blocks have different sizes."//&
01371 CPSourceFileRef,&
01372                error=error)
01373           destination%local_data=REAL(source%local_data_sp,dp)
01374        ELSEIF(.NOT.source%use_sp.AND.destination%use_sp) THEN
01375           CALL cp_assert(SIZE(source%local_data,1)==SIZE(destination%local_data_sp,1).AND.&
01376                SIZE(source%local_data,2)==SIZE(destination%local_data_sp,2),cp_failure_level,&
01377                cp_assertion_failed, routineP,&
01378                "Cannot copy full matrix <"//TRIM(source%name)//&
01379                "> to full matrix <"//TRIM(destination%name)//&
01380                ">. The local_data blocks have different sizes."//&
01381 CPSourceFileRef,&
01382                error=error)
01383           destination%local_data_sp=REAL(source%local_data,sp)
01384        ELSE
01385           CALL cp_assert(SIZE(source%local_data,1)==SIZE(destination%local_data,1).AND.&
01386                SIZE(source%local_data,2)==SIZE(destination%local_data,2),cp_failure_level,&
01387                cp_assertion_failed, routineP,&
01388                "Cannot copy full matrix <"//TRIM(source%name)//&
01389                "> to full matrix <"//TRIM(destination%name)//&
01390                ">. The local_data blocks have different sizes."//&
01391 CPSourceFileRef,&
01392                error=error)
01393           CALL dcopy(SIZE(source%local_data,1)*SIZE(source%local_data,2),&
01394                source%local_data(1,1),1,destination%local_data(1,1),1)
01395        ENDIF
01396     ELSE
01397       CALL cp_assert(source%matrix_struct%nrow_global==&
01398            destination%matrix_struct%nrow_global,&
01399            cp_failure_level,&
01400            cp_assertion_failed, routineP, &
01401            "cannot copy between full matrixes of different sizes"//&
01402 CPSourceFileRef,&
01403            error=error)
01404       CALL cp_assert(source%matrix_struct%ncol_global==&
01405            destination%matrix_struct%ncol_global,&
01406            cp_failure_level,&
01407            cp_assertion_failed, routineP, &
01408            "cannot copy between full matrixes of differen sizes"//&
01409 CPSourceFileRef,&
01410            error=error)
01411 #ifdef __SCALAPACK
01412 
01413       IF(source%use_sp.AND.destination%use_sp) THEN
01414          CALL pscopy(source%matrix_struct%nrow_global*&
01415               source%matrix_struct%ncol_global,&
01416               source%local_data_sp(1,1),1,1,source%matrix_struct%descriptor,1,&
01417               destination%local_data_sp(1,1),1,1,destination%matrix_struct%descriptor,1)
01418       ELSEIF(source%use_sp.AND..NOT.destination%use_sp) THEN
01419          CALL pdcopy(source%matrix_struct%nrow_global*&
01420               source%matrix_struct%ncol_global,&
01421               REAL(source%local_data_sp,dp),1,1,source%matrix_struct%descriptor,1,
01422               destination%local_data(1,1),1,1,destination%matrix_struct%descriptor,1)
01423       ELSEIF(.NOT.source%use_sp.AND.destination%use_sp) THEN
01424          CALL pscopy(source%matrix_struct%nrow_global*&
01425               source%matrix_struct%ncol_global,&
01426               REAL(source%local_data,sp),1,1,source%matrix_struct%descriptor,1,
01427               destination%local_data_sp(1,1),1,1,destination%matrix_struct%descriptor,1)
01428       ELSE
01429          CALL pdcopy(source%matrix_struct%nrow_global*&
01430               source%matrix_struct%ncol_global,&
01431               source%local_data(1,1),1,1,source%matrix_struct%descriptor,1,&
01432               destination%local_data(1,1),1,1,destination%matrix_struct%descriptor,1)
01433       ENDIF
01434 
01435 
01436 #else
01437       CPAssert(.FALSE.,cp_failure_level,routineP,error,failure)
01438 #endif
01439     END IF
01440 
01441     CALL timestop(handle)
01442 
01443   END SUBROUTINE cp_fm_to_fm_matrix
01444 
01445 ! *****************************************************************************
01448   SUBROUTINE cp_fm_to_fm_columns(msource,mtarget,ncol,source_start,&
01449                                  target_start)
01450 
01451     TYPE(cp_fm_type), POINTER                :: msource, mtarget
01452     INTEGER, INTENT(IN)                      :: ncol
01453     INTEGER, INTENT(IN), OPTIONAL            :: source_start, target_start
01454 
01455     CHARACTER(LEN=*), PARAMETER :: routineN = 'cp_fm_to_fm_columns', 
01456       routineP = moduleN//':'//routineN
01457 
01458     INTEGER                                  :: handle, i, n, ss, ts
01459     INTEGER, DIMENSION(9)                    :: desca, descb
01460     REAL(KIND=dp), DIMENSION(:, :), POINTER  :: a, b
01461 
01462 ! -------------------------------------------------------------------------
01463 
01464     CALL timeset(routineN,handle)
01465 
01466     ss = 1
01467     ts = 1
01468 
01469     IF (PRESENT(source_start)) ss = source_start
01470     IF (PRESENT(target_start)) ts = target_start
01471 
01472     n = msource%matrix_struct%nrow_global
01473 
01474     a => msource%local_data
01475     b => mtarget%local_data
01476 
01477 #if defined(__SCALAPACK)
01478     desca(:) = msource%matrix_struct%descriptor(:)
01479     descb(:) = mtarget%matrix_struct%descriptor(:)
01480     DO i=0,ncol-1
01481       CALL pdcopy(n,a(1,1),1,ss+i,desca,1,b(1,1),1,ts+i,descb,1)
01482     END DO
01483 #else
01484     CALL dcopy(ncol*n,a(1,ss),1,b(1,ts),1)
01485 #endif
01486 
01487     CALL timestop(handle)
01488 
01489   END SUBROUTINE cp_fm_to_fm_columns
01490 
01491   SUBROUTINE cp_fm_to_fm_submat(msource, mtarget, nrow, ncol, s_firstrow, s_firstcol, t_firstrow, t_firstcol, error)
01492 
01493     TYPE(cp_fm_type), POINTER                :: msource, mtarget
01494     INTEGER, INTENT(IN)                      :: nrow, ncol, s_firstrow, 
01495                                                 s_firstcol, t_firstrow, 
01496                                                 t_firstcol
01497     TYPE(cp_error_type), INTENT(inout)       :: error
01498 
01499     CHARACTER(LEN=*), PARAMETER :: routineN = 'cp_fm_to_fm_submat', 
01500       routineP = moduleN//':'//routineN
01501 
01502     INTEGER                                  :: handle, i, na, nb, ss, ts
01503     INTEGER, DIMENSION(9)                    :: desca, descb
01504     LOGICAL                                  :: failure
01505     REAL(KIND=dp), DIMENSION(:, :), POINTER  :: a, b
01506     TYPE(cp_blacs_env_type), POINTER         :: context
01507 
01508 ! -------------------------------------------------------------------------
01509 
01510     CALL timeset(routineN,handle)
01511 
01512     failure = .FALSE.
01513 
01514     context => msource%matrix_struct%context
01515 
01516     a => msource%local_data
01517     b => mtarget%local_data
01518 
01519 
01520     na = msource%matrix_struct%nrow_global
01521     nb = mtarget%matrix_struct%nrow_global
01522 !    nrow must be <= na and nb
01523     CALL cp_assert(nrow<=na, cp_failure_level,&
01524            cp_assertion_failed, routineP, &
01525            "cannot copy because nrow > number of rows of source matrix"//&
01526 CPSourceFileRef,&
01527            error=error)
01528     CALL cp_assert(nrow<=nb, cp_failure_level,&
01529            cp_assertion_failed, routineP, &
01530            "cannot copy because nrow > number of rows of target matrix"//&
01531 CPSourceFileRef,&
01532            error=error)
01533     na = msource%matrix_struct%ncol_global
01534     nb = mtarget%matrix_struct%ncol_global
01535 !    ncol must be <= na_col and nb_col
01536     CALL cp_assert(ncol<=na, cp_failure_level,&
01537            cp_assertion_failed, routineP, &
01538            "cannot copy because nrow > number of rows of source matrix"//&
01539 CPSourceFileRef,&
01540            error=error)
01541     CALL cp_assert(ncol<=nb, cp_failure_level,&
01542            cp_assertion_failed, routineP, &
01543            "cannot copy because nrow > number of rows of target matrix"//&
01544 CPSourceFileRef,&
01545            error=error)
01546 
01547 #if defined(__SCALAPACK)
01548     desca(:) = msource%matrix_struct%descriptor(:)
01549     descb(:) = mtarget%matrix_struct%descriptor(:)
01550     DO i=0,ncol-1
01551       ss = s_firstcol+i
01552       ts = t_firstcol+i
01553       CALL pdcopy(nrow,a(1,1),s_firstrow,ss,desca,1,b(1,1),t_firstrow,ts,descb,1)
01554     END DO
01555 #else
01556     DO i = 0,ncol-1
01557       ss = s_firstcol+i
01558       ts = t_firstcol+i
01559       CALL dcopy(nrow,a(s_firstrow,ss),1,b(t_firstrow,ts),1)
01560     END DO
01561 #endif
01562 
01563 
01564     CALL timestop(handle)
01565   END SUBROUTINE cp_fm_to_fm_submat
01566 
01567 ! *****************************************************************************
01568   SUBROUTINE cp_fm_add_to_element(matrix,irow_global,icol_global,alpha,error)
01569 
01570     ! Add alpha to the matrix element specified by the global indices
01571     ! irow_global and icol_global
01572 
01573     ! - Creation (05.05.06,MK)
01574 
01575     TYPE(cp_fm_type), POINTER                :: matrix
01576     INTEGER, INTENT(IN)                      :: irow_global, icol_global
01577     REAL(KIND=dp), INTENT(IN)                :: alpha
01578     TYPE(cp_error_type), INTENT(INOUT)       :: error
01579 
01580     CHARACTER(LEN=*), PARAMETER :: routineN = 'cp_fm_add_to_element', 
01581       routineP = moduleN//':'//routineN
01582 
01583     INTEGER                                  :: icol_local, ipcol, iprow, 
01584                                                 irow_local, mypcol, myprow, 
01585                                                 npcol, nprow
01586     INTEGER, DIMENSION(9)                    :: desca
01587     LOGICAL                                  :: failure
01588     REAL(KIND=dp), DIMENSION(:, :), POINTER  :: a
01589     TYPE(cp_blacs_env_type), POINTER         :: context
01590 
01591 ! -------------------------------------------------------------------------
01592 
01593     failure = .FALSE.
01594 
01595     CPPrecondition(ASSOCIATED(matrix),cp_failure_level,routineP,error,failure)
01596 
01597     context => matrix%matrix_struct%context
01598 
01599     myprow = context%mepos(1)
01600     mypcol = context%mepos(2)
01601 
01602     nprow = context%num_pe(1)
01603     npcol = context%num_pe(2)
01604 
01605     a => matrix%local_data
01606 
01607 #if defined(__SCALAPACK)
01608 
01609     desca(:) = matrix%matrix_struct%descriptor(:)
01610 
01611     CALL infog2l(irow_global,icol_global,desca,nprow,npcol,myprow,mypcol,&
01612                  irow_local,icol_local,iprow,ipcol)
01613 
01614     IF ((iprow == myprow).AND.(ipcol == mypcol)) THEN
01615       a(irow_local,icol_local) = a(irow_local,icol_local) + alpha
01616     END IF
01617 
01618 #else
01619 
01620     a(irow_global,icol_global) = a(irow_global,icol_global) + alpha
01621 
01622 #endif
01623 
01624   END SUBROUTINE cp_fm_add_to_element
01625 
01626  SUBROUTINE cp_fm_write_unformatted(fm,unit,error)
01627     TYPE(cp_fm_type), POINTER                :: fm
01628     INTEGER                                  :: unit
01629     TYPE(cp_error_type), INTENT(inout)       :: error
01630 
01631     CHARACTER(len=*), PARAMETER :: routineN = 'cp_fm_write_unformatted', 
01632       routineP = moduleN//':'//routineN
01633 
01634     INTEGER :: handle, i, i_block, icol_local, ictxt_loc, in, info, ipcol, 
01635       iprow, irow_local, istat, j, max_block, mepos, mypcol, myprow, 
01636       ncol_global, npcol, nprow, nrow_global, num_pe, rb,tag
01637     INTEGER, DIMENSION(9)                    :: desc
01638     LOGICAL                                  :: failure
01639     REAL(KIND=dp), DIMENSION(:), POINTER     :: vecbuf
01640     REAL(KIND=dp), DIMENSION(:, :), POINTER  :: newdat
01641     TYPE(cp_para_env_type), POINTER          :: para_env
01642 
01643 #if defined(__SCALAPACK)
01644     INTEGER, EXTERNAL :: numroc
01645 #endif
01646 
01647    failure=.FALSE.
01648    CALL timeset(routineN,handle)
01649    CALL cp_fm_get_info(fm, nrow_global=nrow_global,ncol_global=ncol_global,ncol_block=max_block,&
01650         para_env=para_env,error=error)
01651 
01652 #if defined(__SCALAPACK)
01653     num_pe=para_env%num_pe
01654     mepos =para_env%mepos
01655     rb=nrow_global
01656     tag =0
01657     ! get a new context
01658     ictxt_loc=para_env%group
01659     CALL cp_blacs_gridinit(ictxt_loc,'R',1,num_pe)
01660     CALL cp_blacs_gridinfo(ictxt_loc,nprow,npcol,myprow,mypcol)
01661     CALL descinit(desc,nrow_global,ncol_global,rb,max_block,0,0,ictxt_loc,nrow_global,info)
01662     CPPostcondition(info==0,cp_failure_level,routineP,error,failure)
01663     in=numroc(ncol_global,max_block,mypcol,0,npcol)
01664 
01665     ALLOCATE(newdat(nrow_global,MAX(1,in)))
01666 
01667             ! do the actual scalapack to cols reordering
01668     CALL pdgemr2d(nrow_global,ncol_global,fm%local_data(1,1),1,1,&
01669          fm%matrix_struct%descriptor, &
01670          newdat(1,1),1,1,desc,ictxt_loc)
01671 
01672 
01673    ALLOCATE(vecbuf(nrow_global*max_block),STAT=istat)
01674    CPPostcondition(istat==0,cp_failure_level,routineP,error,failure)
01675    DO i=1,ncol_global,MAX(max_block,1)
01676       i_block=MIN(max_block,ncol_global-i+1)
01677       CALL infog2l(1,i,desc,nprow,npcol,myprow,mypcol,&
01678                irow_local,icol_local,iprow,ipcol)
01679       IF (ipcol==mypcol) THEN
01680          DO j=1,i_block
01681             vecbuf((j-1)*nrow_global+1:nrow_global*j)=newdat(:,icol_local+j-1)
01682          END DO
01683       END IF
01684 
01685       IF (ipcol==0) THEN
01686          ! do nothing
01687       ELSE
01688          IF (ipcol==mypcol) THEN
01689             CALL mp_send(vecbuf(:),0,tag,para_env%group)
01690          ENDIF
01691          IF (mypcol==0) THEN
01692             CALL mp_recv(vecbuf(:),ipcol,tag,para_env%group)
01693          ENDIF
01694       ENDIF
01695 
01696       IF (unit>0)THEN
01697          DO j=1,i_block
01698             WRITE (unit) vecbuf((j-1)*nrow_global+1:nrow_global*j)
01699          END DO
01700       END IF
01701 
01702    END DO
01703    DEALLOCATE(vecbuf,STAT=istat)
01704    CPPostcondition(istat==0,cp_failure_level,routineP,error,failure)
01705 
01706    CALL cp_blacs_gridexit(ictxt_loc)
01707 
01708    DEALLOCATE(newdat,STAT=istat)
01709    CPPostcondition(istat==0,cp_failure_level,routineP,error,failure)
01710 
01711 #else
01712 
01713    IF (unit>0)THEN
01714       DO j=1,ncol_global
01715          WRITE (unit) fm%local_data(:,j)
01716       END DO
01717    END IF
01718 
01719 #endif
01720 
01721    CALL timestop(handle)
01722  END SUBROUTINE cp_fm_write_unformatted
01723 
01724  SUBROUTINE cp_fm_read_unformatted(fm,unit,error)
01725     TYPE(cp_fm_type), POINTER                :: fm
01726     INTEGER                                  :: unit
01727     TYPE(cp_error_type), INTENT(inout)       :: error
01728 
01729     CHARACTER(len=*), PARAMETER :: routineN = 'cp_fm_read_unformatted', 
01730       routineP = moduleN//':'//routineN
01731 
01732     INTEGER                                  :: handle, j, k, max_block, 
01733                                                 n_cols, ncol_global, 
01734                                                 nrow_global, stat
01735     LOGICAL                                  :: failure
01736     REAL(KIND=dp), DIMENSION(:, :), POINTER  :: vecbuf
01737     TYPE(cp_para_env_type), POINTER          :: para_env
01738 
01739    CALL timeset(routineN,handle)
01740    failure=.FALSE.
01741 
01742    CALL cp_fm_get_info(fm, nrow_global=nrow_global,ncol_global=ncol_global,ncol_block=max_block,&
01743         para_env=para_env,error=error)
01744 
01745 
01746 #if defined(__SCALAPACK)
01747 
01748    ! the parallel case could be made more efficient (see cp_fm_write_unformatted)
01749 
01750    ALLOCATE(vecbuf(nrow_global,max_block),stat=stat)
01751    CPPostcondition(stat==0,cp_failure_level,routineP,error,failure)
01752 
01753    DO j=1,ncol_global,max_block
01754 
01755       n_cols=MIN(max_block,ncol_global-j+1)
01756       IF (para_env%mepos==0) THEN
01757          DO k=1,n_cols
01758             READ(unit) vecbuf(:,k)
01759          ENDDO
01760       ENDIF
01761       CALL mp_bcast(vecbuf,0,para_env%group)
01762       CALL cp_fm_set_submatrix(fm,vecbuf,start_row=1,start_col=j,n_cols=n_cols,error=error)
01763 
01764    ENDDO
01765 
01766    DEALLOCATE(vecbuf,stat=stat)
01767    CPPostcondition(stat==0,cp_failure_level,routineP,error,failure)
01768 
01769 #else
01770 
01771    DO j=1,ncol_global
01772       READ (unit) fm%local_data(:,j)
01773    END DO
01774 
01775 #endif
01776 
01777     CALL timestop(handle)
01778 
01779  END SUBROUTINE cp_fm_read_unformatted
01780 
01781 ! *****************************************************************************
01810  FUNCTION  cp_fm_indxg2p( INDXGLOB, NB, IPROC, ISRCPROC, NPROCS ) RESULT(G2P)
01811    INTEGER, INTENT(IN)                      :: INDXGLOB, NB, IPROC, ISRCPROC, NPROCS
01812    INTEGER                                  :: G2P
01813 
01814    CHARACTER(LEN=*), PARAMETER :: routineN = 'cp_fm_indxg2p', 
01815      routineP = moduleN//':'//routineN
01816 
01817 #if defined(__SCALAPACK)
01818    INTEGER, EXTERNAL :: indxg2p
01819 #endif
01820 
01821 #if defined(__SCALAPACK)
01822 
01823    G2P = indxg2p( INDXGLOB, NB, IPROC, ISRCPROC, NPROCS )
01824 
01825 #else
01826 
01827    G2P = 0
01828 
01829 #endif
01830 
01831  END FUNCTION  cp_fm_indxg2p
01832 
01833 ! *****************************************************************************
01861  FUNCTION  cp_fm_indxg2l( INDXGLOB, NB, IPROC, ISRCPROC, NPROCS ) RESULT(G2L)
01862    INTEGER, INTENT(IN)                      :: INDXGLOB, NB, IPROC, ISRCPROC, NPROCS
01863    INTEGER                                  :: G2L
01864 
01865    CHARACTER(LEN=*), PARAMETER :: routineN = 'cp_fm_indxg2l', 
01866      routineP = moduleN//':'//routineN
01867 
01868 #if defined(__SCALAPACK)
01869    INTEGER, EXTERNAL :: indxg2l
01870 #endif
01871 
01872 #if defined(__SCALAPACK)
01873 
01874    G2L = indxg2l( INDXGLOB, NB, IPROC, ISRCPROC, NPROCS )
01875 
01876 #else
01877 
01878    G2L = INDXGLOB
01879 
01880 #endif
01881 
01882  END FUNCTION  cp_fm_indxg2l
01883 
01884 ! *****************************************************************************
01913  FUNCTION  cp_fm_indxl2g( INDXLOC, NB, IPROC, ISRCPROC, NPROCS ) RESULT(L2G)
01914    INTEGER, INTENT(IN)                      :: INDXLOC, NB, IPROC, ISRCPROC, NPROCS 
01915    INTEGER                                  :: L2G
01916 
01917    CHARACTER(LEN=*), PARAMETER :: routineN = 'cp_fm_indxl2g', 
01918      routineP = moduleN//':'//routineN
01919 
01920 #if defined(__SCALAPACK)
01921    INTEGER, EXTERNAL :: indxl2g
01922 #endif
01923 
01924 #if defined(__SCALAPACK)
01925 
01926    L2G = indxl2g( INDXLOC, NB, IPROC, ISRCPROC, NPROCS )
01927 
01928 #else
01929 
01930    L2G = INDXLOC
01931 
01932 #endif
01933 
01934  END FUNCTION  cp_fm_indxl2g
01935 
01936 END MODULE cp_fm_types