TABLE OF CONTENTS
- ABINIT/m_nctk
- ABINIT/write_eig
- m_nctk/ab_define_var
- m_nctk/bail_if_ncerr
- m_nctk/collect_datar
- m_nctk/create_nc_file
- m_nctk/distrib_datar
- m_nctk/ncfdim_t
- m_nctk/nctk_add_etsf_header
- m_nctk/nctk_def_array_list
- m_nctk/nctk_def_basedims
- m_nctk/nctk_def_dim_list
- m_nctk/nctk_def_dpscalars
- m_nctk/nctk_def_iscalars
- m_nctk/nctk_def_one_array
- m_nctk/nctk_def_one_dim
- m_nctk/nctk_def_scalars_type
- m_nctk/nctk_defnwrite_dpvars
- m_nctk/nctk_defnwrite_ivars
- m_nctk/nctk_defwrite_nonana_raman_terms
- m_nctk/nctk_defwrite_nonana_terms
- m_nctk/nctk_defwrite_raman_terms
- m_nctk/nctk_fort_or_ncfile
- m_nctk/nctk_get_dim
- m_nctk/nctk_idname
- m_nctk/nctk_ncify
- m_nctk/nctk_open_create
- m_nctk/nctk_open_modify
- m_nctk/nctk_open_read
- m_nctk/nctk_read_datar
- m_nctk/nctk_set_atomic_units
- m_nctk/nctk_set_collective
- m_nctk/nctk_set_datamode
- m_nctk/nctk_set_default_for_seq
- m_nctk/nctk_set_defmode
- m_nctk/nctk_string_from_occopt
- m_nctk/nctk_test_mpiio
- m_nctk/nctk_try_fort_or_ncfile
- m_nctk/nctk_write_datar
- m_nctk/nctk_write_dpscalars
- m_nctk/nctk_write_ibz
- m_nctk/nctk_write_iscalars
- m_nctk/nctkarr_t
- m_nctk/nctkerr_t
- m_nctk/nctkvar_t
- m_nctk/str2xtype
- m_nctk/var_from_id
- m_nctk/var_from_name
- m_nctk/write_var_netcdf
ABINIT/m_nctk [ Modules ]
NAME
m_nctk
FUNCTION
Tools and wrappers for NETCDF-IO.
COPYRIGHT
Copyright (C) 2008-2022 ABINIT group (MG) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
TODO
Remove create_nc_file, write_var_netcdf, the output of OUT.nc is dangereous because we can create too many dimensions and get nf90_def_dim - NetCDF library returned: NetCDF: NC_MAX_DIMS exceeded Moreover the multiple calls to redef render the IO very inefficient That part should be rationalized!
SOURCE
23 #if defined HAVE_CONFIG_H 24 #include "config.h" 25 #endif 26 27 #include "abi_common.h" 28 29 MODULE m_nctk 30 31 use defs_basis 32 use m_abicore 33 use m_build_info 34 use m_errors 35 use, intrinsic :: iso_c_binding 36 use m_xmpi 37 #ifdef HAVE_NETCDF 38 use netcdf 39 #endif 40 41 use m_fstrings, only : itoa, sjoin, lstrip, char_count, strcat, endswith, startswith, ltoa 42 use m_io_tools, only : pick_aname, delete_file, file_exists 43 use m_yaml, only : DTSET_IDX 44 45 implicit none 46 47 private
ABINIT/write_eig [ Functions ]
NAME
write_eig
FUNCTION
Write the eigenvalues band by band and k point by k point in a NetCDF file format
INPUTS
filname = Filename of the file where the history will be stored
OUTPUT
(only writing)
SOURCE
2824 subroutine write_eig(eigen,fermie,filename,kptns,mband,nband,nkpt,nsppol,& 2825 & shiftfactor_extfpmd) ! Optional arguments 2826 2827 !Arguments ------------------------------------ 2828 !scalars 2829 character(len=fnlen),intent(in) :: filename 2830 integer,intent(in) :: nkpt,nsppol,mband 2831 real(dp),intent(in) :: fermie 2832 real(dp),optional,intent(in) :: shiftfactor_extfpmd 2833 !arrays 2834 integer,intent(in) :: nband(nkpt*nsppol) 2835 real(dp),intent(in) :: eigen(mband*nkpt*nsppol) 2836 real(dp),intent(in) :: kptns(3,nkpt) 2837 2838 !Local variables------------------------------- 2839 !scalars 2840 integer :: ncerr,ncid,ii, cmode 2841 integer :: xyz_id,nkpt_id,mband_id,nsppol_id 2842 integer :: eig_id,fermie_id,kpt_id,nbk_id,nbk 2843 integer :: shiftfactor_extfpmd_id 2844 integer :: ikpt,isppol,nband_k,band_index 2845 real(dp):: convrt 2846 !arrays 2847 integer :: dimEIG(3),dimKPT(2),dimNBK(2) 2848 integer :: count2(2),start2(2) 2849 integer :: count3(3),start3(3) 2850 integer :: dim0(0) 2851 real(dp):: band(mband) 2852 2853 ! ********************************************************************* 2854 2855 #if defined HAVE_NETCDF 2856 2857 convrt=1.0_dp 2858 2859 !1. Create netCDF file 2860 !ncerr = nf90_create(path=trim(filename),cmode=NF90_CLOBBER, ncid=ncid) 2861 cmode = def_cmode_for_seq_create 2862 ncerr = nf90_create(path=trim(filename), cmode=cmode, ncid=ncid) 2863 NCF_CHECK_MSG(ncerr," create netcdf EIG file") 2864 2865 !2. Define dimensions 2866 ncerr = nf90_def_dim(ncid,"xyz",3,xyz_id) 2867 NCF_CHECK_MSG(ncerr," define dimension xyz") 2868 2869 ncerr = nf90_def_dim(ncid,"mband",mband,mband_id) 2870 NCF_CHECK_MSG(ncerr," define dimension mband") 2871 2872 ncerr = nf90_def_dim(ncid,"nkpt",nkpt,nkpt_id) 2873 NCF_CHECK_MSG(ncerr," define dimension nkpt") 2874 2875 ncerr = nf90_def_dim(ncid,"nsppol",nsppol,nsppol_id) 2876 NCF_CHECK_MSG(ncerr," define dimension nsppol") 2877 2878 !Dimensions for EIGENVALUES 2879 dimEIG = (/ mband_id, nkpt_id, nsppol_id /) 2880 !Dimensions for kpoint positions 2881 dimKPT = (/ xyz_id, nkpt_id /) 2882 !Dimensions for number kpoints per band and spin 2883 !dimNBK = (/ nkpt_id, nsppol_id /) 2884 dimNBK = (/ nkpt_id, nsppol_id /) 2885 2886 !3. Define variables 2887 call ab_define_var(ncid,dim0,fermie_id,NF90_DOUBLE,& 2888 & "fermie","Chemical potential","Hartree") 2889 call ab_define_var(ncid, dimEIG, eig_id, NF90_DOUBLE,& 2890 & "Eigenvalues",& 2891 & "Values of eigenvalues",& 2892 & "Hartree") 2893 call ab_define_var(ncid, dimKPT, kpt_id, NF90_DOUBLE,"Kptns",& 2894 & "Positions of K-points in reciprocal space",& 2895 & "Dimensionless") 2896 call ab_define_var(ncid, dimNBK, nbk_id, NF90_INT,"NBandK",& 2897 & "Number of bands per kpoint and Spin",& 2898 & "Dimensionless") 2899 if(present(shiftfactor_extfpmd)) then 2900 call ab_define_var(ncid,dim0,shiftfactor_extfpmd_id,NF90_DOUBLE,& 2901 & "shiftfactor_extfpmd","Extended FPMD shiftfactor","Hartree") 2902 end if 2903 2904 !4. End define mode 2905 ncerr = nf90_enddef(ncid) 2906 NCF_CHECK_MSG(ncerr," end define mode") 2907 2908 !5 Write kpoint positions 2909 do ikpt=1,nkpt 2910 start2 = (/ 1, ikpt /) 2911 count2 = (/ 3, 1 /) 2912 ncerr = nf90_put_var(ncid, kpt_id,& 2913 & kptns(1:3,ikpt),& 2914 & start = start2,& 2915 & count = count2) 2916 NCF_CHECK_MSG(ncerr," write variable kptns") 2917 end do 2918 2919 !6.1 Write chemical potential 2920 ncerr = nf90_put_var(ncid, fermie_id, fermie) 2921 NCF_CHECK_MSG(ncerr," write variable fermie") 2922 2923 !6.2 Write extfpmd shiftfactor 2924 if(present(shiftfactor_extfpmd)) then 2925 ncerr = nf90_put_var(ncid, shiftfactor_extfpmd_id, shiftfactor_extfpmd) 2926 NCF_CHECK_MSG(ncerr," write variable shiftfactor_extfpmd") 2927 end if 2928 2929 !6.3 Write eigenvalues 2930 band_index=0 2931 do isppol=1,nsppol 2932 do ikpt=1,nkpt 2933 nband_k=nband(ikpt+(isppol-1)*nkpt) 2934 start3 = (/ 1, ikpt, isppol /) 2935 count3 = (/ mband, 1, 1 /) 2936 band(:)=zero 2937 do ii=1,nband_k 2938 band(ii)=eigen(band_index+ii) 2939 end do 2940 ncerr = nf90_put_var(ncid, eig_id,& 2941 & band,& 2942 & start = start3,& 2943 & count = count3) 2944 NCF_CHECK_MSG(ncerr," write variable band") 2945 2946 band_index=band_index+nband_k 2947 end do 2948 end do 2949 2950 !6.4 Write Number of bands per kpoint and Spin 2951 2952 do isppol=1,nsppol 2953 do ikpt=1,nkpt 2954 start2 = (/ ikpt, 1 /) 2955 count2 = (/ 1, 1 /) 2956 nbk=nband(ikpt+(isppol-1)*nkpt) 2957 ncerr = nf90_put_var(ncid, nbk_id,& 2958 & nbk,& 2959 & start = start2) 2960 NCF_CHECK_MSG(ncerr," write variable nband") 2961 end do 2962 end do 2963 2964 !7 Close file 2965 2966 ncerr = nf90_close(ncid) 2967 NCF_CHECK_MSG(ncerr," close netcdf EIG file") 2968 #endif 2969 2970 end subroutine write_eig
m_nctk/ab_define_var [ Functions ]
[ Top ] [ m_nctk ] [ Functions ]
NAME
ab_define_var
FUNCTION
Write the definition of a variable, including units and mnemonics
INPUTS
ncid = Identifier of the netcdf dataset var_dim_id = Identifier of the Dimensions var_id = Identifier of the variable var_mnemo = String of mnemonics var_name = String with the name of the variable var_type = NetCDF type of variable (NF90_DOUBLE, etc) var_units = String of units
OUTPUT
(only writing)
SOURCE
1283 subroutine ab_define_var(ncid, var_dim_id, var_id, var_type, var_name, var_mnemo, var_units) 1284 1285 !Arguments ------------------------------------ 1286 !scalars 1287 integer, intent(in) :: ncid 1288 integer, intent(out) :: var_id 1289 character(len=*), intent(in) :: var_mnemo,var_units,var_name 1290 integer,intent(in) :: var_type 1291 !arrays 1292 integer,intent(in) :: var_dim_id(:) 1293 1294 !Local variables------------------------------- 1295 !scalars 1296 integer :: ncerr 1297 1298 ! ************************************************************************* 1299 1300 #ifdef HAVE_NETCDF 1301 ncerr = nf90_def_var(ncid, trim(var_name), var_type, var_dim_id, var_id) 1302 NCF_CHECK_MSG(ncerr," define variable "//trim(var_name)) 1303 1304 ncerr = nf90_put_att(ncid, var_id, "units",trim(var_units)) 1305 NCF_CHECK_MSG(ncerr," define attribute for "//trim(var_name)) 1306 1307 ncerr = nf90_put_att(ncid, var_id, "mnemonics", trim(var_mnemo)) 1308 NCF_CHECK_MSG(ncerr," define attribute for "//trim(var_name)) 1309 #endif 1310 1311 end subroutine ab_define_var
m_nctk/bail_if_ncerr [ Functions ]
[ Top ] [ m_nctk ] [ Functions ]
NAME
bail_if_ncerr
FUNCTION
INPUTS
ncerr=Netcdf error line=line number of the file where problem occurred file=name of the f90 file containing the caller
SOURCE
642 logical function bail_if_ncerr(ncerr, file, line) result(bail) 643 644 !Arguments ------------------------------------ 645 integer,intent(in) :: ncerr 646 character(len=*),optional,intent(in) :: file 647 integer,optional,intent(in) :: line 648 649 ! ********************************************************************* 650 651 bail = (ncerr /= nf90_noerr) 652 653 if (bail) then 654 einfo%ncerr = ncerr 655 einfo%file = "Subroutine Unknown"; if (present(file)) einfo%file = file 656 einfo%line = 0; if (present(line)) einfo%line = line 657 ! Append Netcdf string to user-defined message. 658 write(einfo%msg,'(2a)')'NetCDF library raised: ',trim(nf90_strerror(ncerr)) 659 end if 660 661 end function bail_if_ncerr
m_nctk/collect_datar [ Functions ]
[ Top ] [ m_nctk ] [ Functions ]
NAME
collect_datar
FUNCTION
Collect a real-space MPI-FFT distributed array on each proc.
INPUTS
ngfft(18)=contain all needed information about 3D FFT (see NOTES at beginning of scfcv) cplex=1 if real array, 2 for complex nfft=Number of FFT points treated by this MPI proc nspden=Second dimension of rhor rhor(cplex*nfft,nspden)=Array in real space (MPI-FFT distributed) fftn3_distrib(n3)=rank of the processors which own fft planes in 3rd dimension. fftn3_local(n3)=local i3 indices comm_fft=MPI-FFT communicator [master]=MPI rank, Optional. If present, the global array is available only on master node.
OUTPUT
rhor_glob(cplex*nfft_tot,nspden)=Global array
SOURCE
2280 subroutine collect_datar(ngfft,cplex,nfft,nspden,rhor,comm_fft,fftn3_distrib,ffti3_local,rhor_glob,master) 2281 2282 !Arguments ------------------------------------ 2283 !scalars 2284 integer,intent(in) :: cplex,nfft,nspden,comm_fft 2285 integer,optional,intent(in) :: master 2286 !arrays 2287 integer,intent(in) :: ngfft(18) 2288 integer,intent(in) :: fftn3_distrib(ngfft(3)),ffti3_local(ngfft(3)) 2289 real(dp),intent(in) :: rhor(cplex*nfft,nspden) 2290 real(dp),intent(out) :: rhor_glob(cplex*product(ngfft(1:3)),nspden) 2291 2292 !Local variables------------------------------- 2293 integer :: ispden,i1,i2,i3,me_fft,i3_local,my_fftbase,glob_fftbase 2294 integer :: n1,n2,n3,ierr,nfft_tot 2295 2296 ! ************************************************************************* 2297 2298 nfft_tot = product(ngfft(1:3)); me_fft = xmpi_comm_rank(comm_fft) 2299 n1 = ngfft(1); n2 = ngfft(2); n3 = ngfft(3) 2300 2301 if (nfft_tot == nfft) then 2302 ! full rhor on each node, just do a copy 2303 rhor_glob = rhor 2304 else 2305 ! if MPI-FFT we have to gather the global array on each node. 2306 rhor_glob = zero 2307 do ispden=1,nspden 2308 do i3=1,n3 2309 if (me_fft == fftn3_distrib(i3)) then 2310 i3_local = ffti3_local(i3) 2311 do i2=1,n2 2312 my_fftbase = cplex * ( (i2-1)*n1 + (i3_local-1)*n1*n2 ) 2313 glob_fftbase = cplex * ( (i2-1)*n1 + (i3-1)*n1*n2 ) 2314 do i1=1,cplex * n1 2315 rhor_glob(i1+glob_fftbase,ispden) = rhor(i1+my_fftbase,ispden) 2316 end do 2317 end do 2318 end if 2319 end do 2320 end do 2321 if (present(master)) then 2322 call xmpi_sum_master(rhor_glob,master,comm_fft,ierr) 2323 else 2324 call xmpi_sum(rhor_glob,comm_fft,ierr) 2325 end if 2326 end if 2327 2328 end subroutine collect_datar
m_nctk/create_nc_file [ Functions ]
[ Top ] [ m_nctk ] [ Functions ]
NAME
create_nc_file
FUNCTION
Create an NetCDF file including a dimension one definition
INPUTS
OUTPUT
TODO: Remove
SOURCE
2692 subroutine create_nc_file (filename,ncid) 2693 2694 !Arguments ------------------------------------ 2695 !scalars 2696 integer,intent(out) :: ncid 2697 !arrays 2698 character(len=*),intent(in) :: filename 2699 2700 !Local variables------------------------------- 2701 #if defined HAVE_NETCDF 2702 integer :: one_id 2703 integer :: ncerr, cmode 2704 #endif 2705 2706 ! ************************************************************************* 2707 2708 ncid = 0 2709 #if defined HAVE_NETCDF 2710 ! Create the NetCDF file 2711 !ncerr = nf90_create(path=filename,cmode=NF90_CLOBBER,ncid=ncid) 2712 cmode = def_cmode_for_seq_create 2713 ncerr = nf90_create(path=filename, cmode=cmode, ncid=ncid) 2714 NCF_CHECK_MSG(ncerr, sjoin('Error while creating:', filename)) 2715 ncerr=nf90_def_dim(ncid,'one',1,one_id) 2716 NCF_CHECK_MSG(ncerr,'nf90_def_dim') 2717 #endif 2718 2719 end subroutine create_nc_file
m_nctk/distrib_datar [ Functions ]
[ Top ] [ m_nctk ] [ Functions ]
NAME
distrib_datar
FUNCTION
distribute a real-space MPI-FFT
INPUTS
ngfft(18)=contain all needed information about 3D FFT (see NOTES at beginning of scfcv) cplex=1 if real array, 2 for complex nfft=Number of FFT points treated by this MPI proc nspden=Second dimension of rhor rhor_glob(cplex*nfft_tot,nspden)=Global array master=The rank of the node that owns the global array. comm_fft=MPI-FFT communicator fftn3_distrib(n3)=rank of the processors which own fft planes in 3rd dimension. fftn3_local(n3)=local i3 indices
OUTPUT
rhor(cplex*nfft,nspden)=Array in real space (MPI-FFT distributed)
SOURCE
2356 subroutine distrib_datar(ngfft,cplex,nfft,nspden,rhor_glob,master,comm_fft,fftn3_distrib,ffti3_local,rhor) 2357 2358 !Arguments ------------------------------------ 2359 !scalars 2360 integer,intent(in) :: cplex,nfft,nspden,comm_fft,master 2361 !arrays 2362 integer,intent(in) :: ngfft(18) 2363 integer,intent(in) :: fftn3_distrib(ngfft(3)),ffti3_local(ngfft(3)) 2364 real(dp),intent(out) :: rhor(cplex*nfft,nspden) 2365 real(dp),intent(inout) :: rhor_glob(cplex*product(ngfft(1:3)),nspden) 2366 2367 !Local variables------------------------------- 2368 integer :: ispden,i1,i2,i3,me_fft,i3_local,my_fftbase,glob_fftbase 2369 integer :: n1,n2,n3,ierr,nfft_tot 2370 2371 ! ************************************************************************* 2372 2373 nfft_tot = product(ngfft(1:3)); me_fft = xmpi_comm_rank(comm_fft) 2374 n1 = ngfft(1); n2 = ngfft(2); n3 = ngfft(3) 2375 2376 if (nfft_tot == nfft) then 2377 ! full rhor on each node, just do a copy 2378 rhor = rhor_glob 2379 else 2380 ! if MPI-FFT we have to gather the global array on each node. 2381 call xmpi_bcast(rhor_glob,master,comm_fft,ierr) 2382 do ispden=1,nspden 2383 do i3=1,n3 2384 if (me_fft == fftn3_distrib(i3)) then 2385 i3_local = ffti3_local(i3) 2386 do i2=1,n2 2387 my_fftbase = cplex * ( (i2-1)*n1 + (i3_local-1)*n1*n2 ) 2388 glob_fftbase = cplex * ( (i2-1)*n1 + (i3-1)*n1*n2 ) 2389 do i1=1,cplex * n1 2390 rhor(i1+my_fftbase,ispden) = rhor_glob(i1+glob_fftbase,ispden) 2391 end do 2392 end do 2393 end if 2394 end do 2395 end do 2396 end if 2397 2398 end subroutine distrib_datar
m_nctk/ncfdim_t [ Types ]
NAME
nctkdim_t
FUNCTION
Stores the name and the value of a netcdf dimension
SOURCE
119 type,public :: nctkdim_t 120 character(len=nctk_slen) :: name ! name of the dimension. 121 integer :: value ! value of the dimension. 122 !integer :: id 123 end type nctkdim_t
m_nctk/nctk_add_etsf_header [ Functions ]
[ Top ] [ m_nctk ] [ Functions ]
NAME
nctk_add_etsf_header
FUNCTION
Add the etsf-io header to a file associated to a netcdf file handler.
INPUTS
ncid=Netcdf file identifier. * version = the number of version to be created. * title = (optional) a title for the file (80 characters max). * history = (optional) the first line of history (1024 characters max). * with_etsf_header = (optional) if true, will create a header as defined in the ETSF specifications (default is .true.). When value is .false., arguments title, history and version are ignored.
SOURCE
878 integer function nctk_add_etsf_header(ncid, title, history) result(ncerr) 879 880 !Arguments ------------------------------------ 881 integer,intent(in) :: ncid 882 character(len=*),optional,intent(in) :: title,history 883 884 !Local variables------------------------------- 885 !integer :: ncerr 886 character(len=*), parameter :: file_format = "ETSF Nanoquanta" 887 character(len=*), parameter :: conventions = "http://www.etsf.eu/fileformats/" 888 real :: format_version = 3.3 ! Real is not a good choice for a version! 889 890 ! ********************************************************************* 891 892 ncerr = nctk_set_defmode(ncid) 893 if (ncerr /= nf90_noerr) return 894 895 ! The file format 896 ncerr = nf90_put_att(ncid, NF90_GLOBAL, "file_format", file_format) 897 if (ncerr /= nf90_noerr) return 898 899 ! The version 900 ncerr = nf90_put_att(ncid, NF90_GLOBAL, "file_format_version", format_version) 901 if (ncerr /= nf90_noerr) return 902 903 ! The conventions 904 ncerr = nf90_put_att(ncid, NF90_GLOBAL, "Conventions", conventions) 905 if (ncerr /= nf90_noerr) return 906 907 ! The history if present 908 if (present(history)) then 909 ncerr = nf90_put_att(ncid, NF90_GLOBAL, "history", history(1:min(1024, len(history)))) 910 if (ncerr /= nf90_noerr) return 911 end if 912 913 ! The title if present 914 if (present(title)) then 915 ncerr = nf90_put_att(ncid, NF90_GLOBAL, "title", title(1:min(80, len(title)))) 916 if (ncerr /= nf90_noerr) return 917 end if 918 919 ! These are extensions not in the standard. 920 ! Add info on the code that produced this file 921 ncerr = nf90_put_att(ncid, NF90_GLOBAL, "code", "Abinit") 922 if (ncerr /= nf90_noerr) return 923 924 ncerr = nf90_put_att(ncid, NF90_GLOBAL, "code_version", ABINIT_VERSION) 925 if (ncerr /= nf90_noerr) return 926 927 end function nctk_add_etsf_header
m_nctk/nctk_def_array_list [ Functions ]
[ Top ] [ m_nctk ] [ Functions ]
NAME
nctk_def_array_list
FUNCTION
Define list of arrays with a given type, Return immediately if error
INPUTS
ncid=Netcdf identifier. nctk_arrays(:)=List of array descriptors. [defmode]=If True, the nc file is set in define mode (default=False) [prefix]=Prefix added to varnames and dimensions. Empty string if not specified.
SOURCE
1599 integer function nctk_def_array_list(ncid, nctk_arrays, defmode, prefix) result(ncerr) 1600 1601 !Arguments ------------------------------------ 1602 !scalars 1603 integer,intent(in) :: ncid 1604 logical,optional,intent(in) :: defmode 1605 character(len=*),optional,intent(in) :: prefix 1606 !arrays 1607 type(nctkarr_t),intent(in) :: nctk_arrays(:) 1608 1609 !Local variables------------------------------- 1610 !scalars 1611 integer :: ia 1612 1613 ! ********************************************************************* 1614 1615 ncerr = nf90_noerr 1616 if (present(defmode)) then 1617 if (defmode) then 1618 NCF_CHECK(nctk_set_defmode(ncid)) 1619 end if 1620 end if 1621 1622 do ia=1,size(nctk_arrays) 1623 if (present(prefix)) then 1624 NCF_CHECK(nctk_def_one_array(ncid, nctk_arrays(ia), prefix=prefix)) 1625 else 1626 NCF_CHECK(nctk_def_one_array(ncid, nctk_arrays(ia))) 1627 end if 1628 end do 1629 1630 end function nctk_def_array_list
m_nctk/nctk_def_basedims [ Functions ]
[ Top ] [ m_nctk ] [ Functions ]
NAME
nctk_def_basedims
FUNCTION
Define the basic dimensions used in ETSF-IO files.
INPUTS
ncid=Netcdf identifier. [defmode]=If True, the nc file is set in define mode (default=False)
SOURCE
1225 integer function nctk_def_basedims(ncid, defmode) result(ncerr) 1226 1227 !Arguments ------------------------------------ 1228 integer,intent(in) :: ncid 1229 logical,optional,intent(in) :: defmode 1230 1231 ! ********************************************************************* 1232 1233 ncerr = nf90_noerr 1234 1235 if (present(defmode)) then 1236 if (defmode) then 1237 NCF_CHECK(nctk_set_defmode(ncid)) 1238 end if 1239 end if 1240 1241 ! Basic ETSF-IO dimensions that should be always present in the file. 1242 ncerr = nctk_def_dims(ncid, [& 1243 nctkdim_t("complex", 2), nctkdim_t("symbol_length", 2), nctkdim_t("character_string_length", etsfio_charlen),& 1244 nctkdim_t("number_of_cartesian_directions", 3), nctkdim_t("number_of_reduced_dimensions", 3),& 1245 nctkdim_t("number_of_vectors", 3) & 1246 ]) 1247 NCF_CHECK(ncerr) 1248 1249 ! Useful integers. 1250 ncerr = nctk_def_dims(ncid, [ & 1251 nctkdim_t("one", 1), nctkdim_t("two", 2), nctkdim_t("three", 3), & 1252 nctkdim_t("four", 4), nctkdim_t("five", 5), nctkdim_t("six", 6), & 1253 nctkdim_t("seven", 7), nctkdim_t("eight", 8), nctkdim_t("nine", 9), nctkdim_t("ten", 10), & 1254 nctkdim_t("fnlen", fnlen + 1) & 1255 ]) 1256 NCF_CHECK(ncerr) 1257 1258 end function nctk_def_basedims
m_nctk/nctk_def_dim_list [ Functions ]
[ Top ] [ m_nctk ] [ Functions ]
NAME
nctk_def_dim_list
FUNCTION
Define list of dimensions variables and write their values. Return immediately if error
INPUTS
ncid=Netcdf identifier. dimnames(:)=List of strings with the name of the dimensions values(:)=List of integer scalars [defmode]=If True, the nc file is set in define mode (default=False) [prefix]=Prefix added to varnames and dimensions. Empty string if not specified.
SOURCE
1139 integer function nctk_def_dim_list(ncid, nctkdims, defmode, prefix) result(ncerr) 1140 1141 !Arguments ------------------------------------ 1142 !scalars 1143 integer,intent(in) :: ncid 1144 logical,optional,intent(in) :: defmode 1145 character(len=*),optional,intent(in) :: prefix 1146 !arrays 1147 type(nctkdim_t),intent(in) :: nctkdims(:) 1148 1149 !Local variables------------------------------- 1150 !scalars 1151 integer :: ii 1152 1153 ! ********************************************************************* 1154 1155 ncerr = nf90_noerr 1156 if (present(defmode)) then 1157 if (defmode) then 1158 NCF_CHECK(nctk_set_defmode(ncid)) 1159 end if 1160 end if 1161 1162 do ii=1,size(nctkdims) 1163 if (present(prefix)) then 1164 ncerr = nctk_def_one_dim(ncid, nctkdims(ii), prefix=prefix) 1165 else 1166 ncerr = nctk_def_one_dim(ncid, nctkdims(ii)) 1167 end if 1168 if (ncerr /= nf90_noerr) return 1169 end do 1170 1171 end function nctk_def_dim_list
m_nctk/nctk_def_dpscalars [ Functions ]
[ Top ] [ m_nctk ] [ Functions ]
NAME
nctk_def_dpscalars
FUNCTION
Define list of double precision **scalar** variables. Return immediately if error
INPUTS
ncid=Netcdf identifier. varnames(:)=List of strings with the name of the variables [defmode]=If True, the nc file is set in define mode (default=False) [prefix]=Prefix added to varnames and dimensions. Empty string if not specified.
SOURCE
1435 integer function nctk_def_dpscalars(ncid, varnames, defmode, prefix) result(ncerr) 1436 1437 !Arguments ------------------------------------ 1438 !scalars 1439 integer,intent(in) :: ncid 1440 logical,optional,intent(in) :: defmode 1441 character(len=*),optional,intent(in) :: prefix 1442 !arrays 1443 character(len=*),intent(in) :: varnames(:) 1444 1445 !Local variables------------------------------- 1446 character(len=nctk_slen) :: prefix_ 1447 1448 ! ********************************************************************* 1449 prefix_ = ""; if (present(prefix)) prefix_ = prefix 1450 1451 if (present(defmode)) then 1452 ncerr = nctk_def_scalars_type(ncid, varnames, nf90_double, defmode=defmode, prefix=prefix_) 1453 else 1454 ncerr = nctk_def_scalars_type(ncid, varnames, nf90_double, prefix=prefix_) 1455 end if 1456 1457 end function nctk_def_dpscalars
m_nctk/nctk_def_iscalars [ Functions ]
[ Top ] [ m_nctk ] [ Functions ]
NAME
nctk_def_iscalars
FUNCTION
Define list of integer **scalar** variables. Return immediately if error
INPUTS
ncid=Netcdf identifier. varnames(:)=List of strings with the name of the variables [defmode]=If True, the nc file is set in define mode (default=False) [prefix]=Prefix added to varnames and dimensions. Empty string if not specified.
SOURCE
1395 integer function nctk_def_iscalars(ncid, varnames, defmode, prefix) result(ncerr) 1396 1397 !Arguments ------------------------------------ 1398 !scalars 1399 integer,intent(in) :: ncid 1400 logical,optional,intent(in) :: defmode 1401 character(len=*),optional,intent(in) :: prefix 1402 !arrays 1403 character(len=*),intent(in) :: varnames(:) 1404 1405 ! ********************************************************************* 1406 1407 if (present(defmode)) then 1408 ncerr = nctk_def_scalars_type(ncid, varnames, nf90_int, defmode=defmode) 1409 else 1410 if (present(prefix)) then 1411 ncerr = nctk_def_scalars_type(ncid, varnames, nf90_int, prefix=prefix) 1412 else 1413 ncerr = nctk_def_scalars_type(ncid, varnames, nf90_int) 1414 end if 1415 end if 1416 1417 end function nctk_def_iscalars
m_nctk/nctk_def_one_array [ Functions ]
[ Top ] [ m_nctk ] [ Functions ]
NAME
nctk_def_one_array
FUNCTION
Define list of arrays with a given type, Return immediately if error
INPUTS
ncid=Netcdf identifier. nctk_array=Array descriptor. [defmode]=If True, the nc file is set in define mode (default=False) [prefix]=Prefix added to varnames and dimensions. Empty string if not specified.
SOURCE
1476 integer function nctk_def_one_array(ncid, nctk_array, defmode, varid, prefix) result(ncerr) 1477 1478 !Arguments ------------------------------------ 1479 !scalars 1480 integer,intent(in) :: ncid 1481 logical,optional,intent(in) :: defmode 1482 integer,optional,intent(out) :: varid 1483 type(nctkarr_t),intent(in) :: nctk_array 1484 character(len=*),optional,intent(in) :: prefix 1485 1486 !Local variables------------------------------- 1487 !scalars 1488 integer :: ii,xtype,prev,cnt,nn,vid 1489 character(len=500) :: msg 1490 !arrays 1491 integer :: dimids(NF90_MAX_DIMS),dimvals(NF90_MAX_DIMS) 1492 character(len=nctk_slen) :: sarr(NF90_MAX_DIMS), string, pre, vname, dimname 1493 type(nctkvar_t) :: var 1494 1495 ! ********************************************************************* 1496 1497 pre = ""; if (present(prefix)) pre = prefix 1498 1499 ncerr = nf90_noerr 1500 if (present(defmode)) then 1501 if (defmode) then 1502 NCF_CHECK(nctk_set_defmode(ncid)) 1503 end if 1504 end if 1505 1506 xtype = str2xtype(nctk_array%dtype) 1507 vname = strcat(pre, nctk_array%name) 1508 1509 ! Build array of strings with the dimensions. 1510 string = nctk_array%shape_str 1511 nn = char_count(string, ",") 1512 ABI_CHECK(nn <= NF90_MAX_DIMS, "Too many dimensions!") 1513 1514 ! Parse dimension names and add prefix (if any). 1515 if (nn == 0) then 1516 cnt = 1 1517 dimname = lstrip(string) 1518 if (any(dimname == NCTK_IMPLICIT_DIMS)) then 1519 sarr(1) = dimname 1520 else 1521 sarr(1) = strcat(pre, dimname) 1522 end if 1523 else 1524 prev = 0; cnt = 0 1525 do ii=1,len_trim(string) 1526 if (string(ii:ii) == ",") then 1527 cnt = cnt + 1 1528 dimname = lstrip(string(prev+1:ii-1)) 1529 if (any(dimname == NCTK_IMPLICIT_DIMS)) then 1530 sarr(cnt) = dimname 1531 else 1532 sarr(cnt) = strcat(pre, dimname) 1533 end if 1534 prev = ii 1535 end if 1536 end do 1537 cnt = cnt + 1 1538 dimname = lstrip(string(prev+1:ii-1)) 1539 if (any(dimname == NCTK_IMPLICIT_DIMS)) then 1540 sarr(cnt) = dimname 1541 else 1542 sarr(cnt) = strcat(pre, dimname) 1543 end if 1544 end if 1545 1546 ! Get dimids 1547 do ii=1,cnt 1548 NCF_CHECK_MSG(nf90_inq_dimid(ncid, sarr(ii), dimids(ii)), sarr(ii)) 1549 NCF_CHECK(nf90_inquire_dimension(ncid, dimids(ii), len=dimvals(ii))) 1550 end do 1551 1552 ! Check if dimension already exists. 1553 ! Variable already exists. Check if type and dimensions agree 1554 if (nf90_inq_varid(ncid, vname, vid) == nf90_noerr) then 1555 call var_from_id(ncid, vid, var) 1556 if (.not. (var%xtype == xtype .and. var%ndims == cnt)) then 1557 write(msg,"(4a,2(2(a,i0),a))")& 1558 "variable ",trim(vname)," already exists with a different definition:",ch10,& 1559 "In file: xtype = ",var%xtype,", ndims = ",var%ndims,ch10,& 1560 "From caller: xtype = ",xtype,", ndims = ",cnt,ch10 1561 ABI_ERROR(msg) 1562 end if 1563 if (any(dimvals(1:cnt) /= var%dimlens(1:var%ndims))) then 1564 write(msg,"(4a,2(3a))")& 1565 "variable ",trim(vname)," already exists but with different shape.",ch10,& 1566 "In file: dims = ",trim(ltoa(var%dimlens(:var%ndims))),ch10,& 1567 "From caller dims = ",trim(ltoa(dimvals(:cnt))),ch10 1568 ABI_ERROR(msg) 1569 end if 1570 if (present(varid)) varid = vid 1571 return 1572 end if 1573 1574 ! Define variable since it doesn't exist. 1575 ncerr = nf90_def_var(ncid, vname, xtype, dimids(1:cnt), vid) 1576 NCF_CHECK(ncerr) 1577 1578 if (present(varid)) varid = vid 1579 1580 end function nctk_def_one_array
m_nctk/nctk_def_one_dim [ Functions ]
[ Top ] [ m_nctk ] [ Functions ]
NAME
nctk_def_one_dim
FUNCTION
Define list of dimensions variables and write their values. Return immediately if error
INPUTS
ncid=Netcdf identifier. dimnames(:)=List of strings with the name of the dimensions values(:)=List of integer scalars [defmode]=If True, the nc file is set in define mode (default=False) [prefix]=Prefix added to varnames and dimensions. Empty string if not specified.
SOURCE
1067 integer function nctk_def_one_dim(ncid, nctkdim, defmode, prefix) result(ncerr) 1068 1069 !Arguments ------------------------------------ 1070 !scalars 1071 integer,intent(in) :: ncid 1072 logical,optional,intent(in) :: defmode 1073 character(len=*),optional,intent(in) :: prefix 1074 !arrays 1075 type(nctkdim_t),intent(in) :: nctkdim 1076 1077 !Local variables------------------------------- 1078 integer :: dimid,dimlen 1079 character(len=nctk_slen) :: dname 1080 character(len=500) :: msg 1081 ! ********************************************************************* 1082 1083 ncerr = nf90_noerr 1084 1085 if (present(defmode)) then 1086 if (defmode) then 1087 NCF_CHECK(nctk_set_defmode(ncid)) 1088 end if 1089 end if 1090 1091 if (present(prefix)) then 1092 if (any(nctkdim%name == NCTK_IMPLICIT_DIMS)) then 1093 dname = nctkdim%name 1094 else 1095 dname = strcat(prefix, nctkdim%name) 1096 end if 1097 else 1098 dname = nctkdim%name 1099 end if 1100 1101 ! if dimension already exists, test whether it has the same value else define new dim. 1102 ncerr = nf90_inq_dimid(ncid, dname, dimid) 1103 1104 if (ncerr == nf90_noerr) then 1105 NCF_CHECK(nf90_inquire_dimension(ncid, dimid, len=dimlen)) 1106 if (dimlen /= nctkdim%value) then 1107 write(msg, "(4a,2(a,i0))")& 1108 "dimension ", trim(dname)," already exists but with a different value",ch10,& 1109 "from file: ", dimlen, "; about to write: ", nctkdim%value 1110 ABI_ERROR(msg) 1111 end if 1112 else 1113 ncerr = nf90_def_dim(ncid, dname, nctkdim%value, dimid) 1114 NCF_CHECK(ncerr) 1115 end if 1116 1117 end function nctk_def_one_dim
m_nctk/nctk_def_scalars_type [ Functions ]
[ Top ] [ m_nctk ] [ Functions ]
NAME
nctk_def_scalars_type
FUNCTION
Define list of **scalar** variables with a given type, Return immediately if error
INPUTS
ncid=Netcdf identifier. varnames(:)=List of strings with the name of the variables xtype=Type of the variables [defmode]=If True, the nc file is set in define mode (default=False) [prefix]=Prefix added to varnames and dimensions. Empty string if not specified.
SOURCE
1330 integer function nctk_def_scalars_type(ncid, varnames, xtype, defmode, prefix) result(ncerr) 1331 1332 !Arguments ------------------------------------ 1333 !scalars 1334 integer,intent(in) :: ncid,xtype 1335 logical,optional,intent(in) :: defmode 1336 character(len=*),optional,intent(in) :: prefix 1337 !arrays 1338 character(len=*),intent(in) :: varnames(:) 1339 1340 !Local variables------------------------------- 1341 !scalars 1342 integer :: ii,varid 1343 character(len=nctk_slen) :: vname 1344 type(nctkvar_t) :: var 1345 1346 ! ********************************************************************* 1347 1348 ncerr = nf90_noerr 1349 if (present(defmode)) then 1350 if (defmode) then 1351 NCF_CHECK(nctk_set_defmode(ncid)) 1352 end if 1353 end if 1354 1355 ! Special case where dimension is null 1356 do ii=1,size(varnames) 1357 vname = varnames(ii) 1358 if (present(prefix)) vname = strcat(prefix, vname) 1359 1360 if (nf90_inq_varid(ncid, vname, varid) == nf90_noerr) then 1361 ! Variable already exists. Check if type and dimensions agree 1362 call var_from_id(ncid, varid, var) 1363 1364 if (.not. (var%xtype == xtype .and. var%ndims == 0)) then 1365 ABI_ERROR("variable already exists with a different definition.") 1366 else 1367 cycle ! Dimension matches, skip definition. 1368 end if 1369 1370 else 1371 ! Define variable since it doesn't exist. 1372 ncerr = nf90_def_var(ncid, vname, xtype, varid) 1373 NCF_CHECK(ncerr) 1374 end if 1375 end do 1376 1377 end function nctk_def_scalars_type
m_nctk/nctk_defnwrite_dpvars [ Functions ]
[ Top ] [ m_nctk ] [ Functions ]
NAME
nctk_defnwrite_dpvars
FUNCTION
Define list of real(dp) **scalar** variables and write their values. Return immediately if error
INPUTS
ncid=Netcdf identifier. varnames(:)=List of strings with the name of the variables values(:)=List of integer scalars
SOURCE
1792 integer function nctk_defnwrite_dpvars(ncid, varnames, values) result(ncerr) 1793 1794 !Arguments ------------------------------------ 1795 !scalars 1796 integer,intent(in) :: ncid 1797 !arrays 1798 real(dp),intent(in) :: values(:) 1799 character(len=*),intent(in) :: varnames(:) 1800 1801 !Local variables------------------------------- 1802 !scalars 1803 integer :: ii,varid 1804 !arrays 1805 1806 ! ********************************************************************* 1807 ncerr = nf90_noerr 1808 1809 ABI_CHECK(size(varnames) == size(values), "Different size in varnames, values") 1810 1811 ncerr = nctk_def_dpscalars(ncid, varnames, defmode=.True.) 1812 NCF_CHECK(ncerr) 1813 1814 NCF_CHECK(nctk_set_datamode(ncid)) 1815 do ii=1,size(varnames) 1816 !write(std_out,*)varnames(ii) 1817 varid = nctk_idname(ncid, varnames(ii)) 1818 NCF_CHECK(nf90_put_var(ncid, varid, values(ii))) 1819 end do 1820 1821 end function nctk_defnwrite_dpvars
m_nctk/nctk_defnwrite_ivars [ Functions ]
[ Top ] [ m_nctk ] [ Functions ]
NAME
nctk_defnwrite_ivars
FUNCTION
Define list of integer **scalar** variables and write their values. Return immediately if error
INPUTS
ncid=Netcdf identifier. varnames(:)=List of strings with the name of the variables values(:)=List of integer scalars
SOURCE
1748 integer function nctk_defnwrite_ivars(ncid, varnames, values) result(ncerr) 1749 1750 !Arguments ------------------------------------ 1751 !scalars 1752 integer,intent(in) :: ncid 1753 !arrays 1754 integer,intent(in) :: values(:) 1755 character(len=*),intent(in) :: varnames(:) 1756 1757 !Local variables------------------------------- 1758 !scalars 1759 integer :: ii,varid 1760 1761 ! ********************************************************************* 1762 1763 ABI_CHECK(size(varnames) == size(values), "Different size in varnames, values") 1764 1765 ncerr = nctk_def_iscalars(ncid, varnames, defmode=.True.) 1766 NCF_CHECK(ncerr) 1767 1768 NCF_CHECK(nctk_set_datamode(ncid)) 1769 do ii=1,size(varnames) 1770 varid = nctk_idname(ncid, varnames(ii)) 1771 NCF_CHECK(nf90_put_var(ncid, varid, values(ii))) 1772 end do 1773 1774 end function nctk_defnwrite_ivars
m_nctk/nctk_defwrite_nonana_raman_terms [ Functions ]
[ Top ] [ m_nctk ] [ Functions ]
NAME
nctk_defwrite_nonana_raman_terms
FUNCTION
Write the Raman susceptiblities for q-->0 along different directions in the netcdf file.
INPUTS
ncid=netcdf file id. iphl2=Index of the q-point to be written to file. nph2l=Number of qpoints. rsus(3*natom,3,3)=List of Raman susceptibilities along the direction corresponding to iphl2. natom=Number of atoms
OUTPUT
Only writing.
SOURCE
2584 subroutine nctk_defwrite_nonana_raman_terms(ncid, iphl2, nph2l, natom, rsus, mode) 2585 2586 !Arguments ------------------------------------ 2587 !scalars 2588 integer,intent(in) :: ncid,natom,iphl2,nph2l 2589 character(len=*),intent(in) :: mode 2590 !arrays 2591 real(dp),intent(in) :: rsus(3*natom,3,3) 2592 2593 !Local variables------------------------------- 2594 !scalars 2595 integer :: ncerr, raman_sus_varid 2596 2597 ! ************************************************************************* 2598 2599 !Fake use of nph2l, to keep it as argument. This should be removed when nph2l will be used. 2600 if(.false.)then 2601 ncerr=nph2l 2602 endif 2603 2604 select case (mode) 2605 case ("define") 2606 NCF_CHECK(nctk_def_basedims(ncid, defmode=.True.)) 2607 ncerr = nctk_def_arrays(ncid, [ nctkarr_t("non_analytical_raman_sus", "dp", & 2608 "number_of_non_analytical_directions,number_of_phonon_modes,number_of_cartesian_directions,number_of_cartesian_directions")]) 2609 NCF_CHECK(ncerr) 2610 2611 NCF_CHECK(nctk_set_datamode(ncid)) 2612 2613 case ("write") 2614 NCF_CHECK(nf90_inq_varid(ncid, "non_analytical_raman_sus", raman_sus_varid)) 2615 ncerr = nf90_put_var(ncid,raman_sus_varid,rsus,& 2616 start=[iphl2,1,1,1], count=[1,3*natom,3,3]) 2617 NCF_CHECK(ncerr) 2618 2619 case default 2620 ABI_ERROR(sjoin("Wrong value for mode", mode)) 2621 end select 2622 2623 end subroutine nctk_defwrite_nonana_raman_terms
m_nctk/nctk_defwrite_nonana_terms [ Functions ]
[ Top ] [ m_nctk ] [ Functions ]
NAME
nctk_defwrite_nonana_terms
FUNCTION
Write to ncfile the phonon frequencies and displacements for q --> 0 in the presence of non-analytical behaviour.
INPUTS
ncid=netcdf file id. iphl2=Index of the q-point to be written to file nph2l=Number of qpoints. qph2l(3,nph2l)=List of phonon wavevector directions along which the non-analytical correction to the Gamma-point phonon frequencies will be calculated The direction is in CARTESIAN COORDINATES natom=Number of atoms phfrq(3*natom)=Phonon frequencies in Ha cart_displ(2,3*natom,3*natom)=displacements in CARTESIAN coordinates.
OUTPUT
Only writing.
SOURCE
2516 subroutine nctk_defwrite_nonana_terms(ncid, iphl2, nph2l, qph2l, natom, phfrq, cart_displ, mode) 2517 2518 !Arguments ------------------------------------ 2519 !scalars 2520 integer,intent(in) :: ncid,iphl2,nph2l,natom 2521 character(len=*),intent(in) :: mode 2522 !arrays 2523 real(dp),intent(in) :: qph2l(3, nph2l) 2524 real(dp),intent(in) :: phfrq(3*natom) 2525 real(dp),intent(in) :: cart_displ(2,3*natom,3*natom) 2526 2527 !Local variables------------------------------- 2528 !scalars 2529 integer :: ncerr, na_phmodes_varid, na_phdispl_varid 2530 2531 ! ************************************************************************* 2532 2533 select case (mode) 2534 case ("define") 2535 !NCF_CHECK(nctk_def_basedims(ncid, defmode=.True.)) 2536 ncerr = nctk_def_dims(ncid, [nctkdim_t("number_of_non_analytical_directions", nph2l)], defmode=.True.) 2537 NCF_CHECK(ncerr) 2538 2539 ncerr = nctk_def_arrays(ncid, [& 2540 nctkarr_t('non_analytical_directions', "dp", "number_of_cartesian_directions, number_of_non_analytical_directions"),& 2541 nctkarr_t('non_analytical_phonon_modes', "dp", "number_of_phonon_modes, number_of_non_analytical_directions"),& 2542 nctkarr_t('non_analytical_phdispl_cart', "dp", & 2543 "two, number_of_phonon_modes, number_of_phonon_modes, number_of_non_analytical_directions")]) 2544 NCF_CHECK(ncerr) 2545 2546 NCF_CHECK(nctk_set_datamode(ncid)) 2547 NCF_CHECK(nf90_put_var(ncid, nctk_idname(ncid, "non_analytical_directions"), qph2l)) 2548 2549 case ("write") 2550 2551 NCF_CHECK(nf90_inq_varid(ncid, "non_analytical_phonon_modes", na_phmodes_varid)) 2552 NCF_CHECK(nf90_put_var(ncid,na_phmodes_varid,phfrq*Ha_eV,start=[1, iphl2], count=[3*natom, 1])) 2553 NCF_CHECK(nf90_inq_varid(ncid, "non_analytical_phdispl_cart", na_phdispl_varid)) 2554 ncerr = nf90_put_var(ncid,na_phdispl_varid,cart_displ*Bohr_Ang,& 2555 start=[1,1,1,iphl2], count=[2,3*natom,3*natom, 1]) 2556 NCF_CHECK(ncerr) 2557 2558 case default 2559 ABI_ERROR(sjoin("Wrong value for mode", mode)) 2560 end select 2561 2562 end subroutine nctk_defwrite_nonana_terms
m_nctk/nctk_defwrite_raman_terms [ Functions ]
[ Top ] [ m_nctk ] [ Functions ]
NAME
nctk_defwrite_raman_terms
FUNCTION
Write the Raman susceptiblities for q=0 and also the phonon frequncies at gamma.
INPUTS
ncid=netcdf file id. rsus(3*natom,3,3)=List of Raman susceptibilities. natom=Number of atoms
OUTPUT
Only writing.
SOURCE
2643 subroutine nctk_defwrite_raman_terms(ncid, natom, rsus, phfrq) 2644 2645 !Arguments ------------------------------------ 2646 !scalars 2647 integer,intent(in) :: ncid,natom 2648 !arrays 2649 real(dp),intent(in) :: rsus(3*natom,3,3) 2650 real(dp),intent(in) :: phfrq(3*natom) 2651 2652 !Local variables------------------------------- 2653 !scalars 2654 integer :: ncerr, raman_sus_varid, phmodes_varid 2655 2656 ! ************************************************************************* 2657 2658 NCF_CHECK(nctk_def_basedims(ncid, defmode=.True.)) 2659 ncerr = nctk_def_arrays(ncid, [ nctkarr_t("raman_sus", "dp", & 2660 "number_of_phonon_modes,number_of_cartesian_directions,number_of_cartesian_directions"), & 2661 nctkarr_t("gamma_phonon_modes", "dp", "number_of_phonon_modes")]) 2662 NCF_CHECK(ncerr) 2663 2664 NCF_CHECK(nctk_set_datamode(ncid)) 2665 2666 NCF_CHECK(nf90_inq_varid(ncid, "raman_sus", raman_sus_varid)) 2667 NCF_CHECK(nf90_put_var(ncid,raman_sus_varid,rsus)) 2668 NCF_CHECK(nf90_inq_varid(ncid, "gamma_phonon_modes", phmodes_varid)) 2669 NCF_CHECK(nf90_put_var(ncid,phmodes_varid,phfrq*Ha_eV)) 2670 2671 end subroutine nctk_defwrite_raman_terms
m_nctk/nctk_fort_or_ncfile [ Functions ]
[ Top ] [ m_nctk ] [ Functions ]
NAME
nctk_fort_or_ncfile
FUNCTION
Return the iomode used to perform IO operations on filename. If filename does not exist, a similar file with extension `.nc` is tried and iomode is set to IO_MODE_ETSF if the file exists. This trick is used to run the Abinit test suite in netcdf mode without changing the input files. The modification (if any) is logged to std_out.
SIDE EFFECTS
filename=Tentative filename in input. Changed to netcdf file if input filename does not exist and a file with netcdf extension is found.
OUTPUT
iomode=Flag selecting the IO library. Set to IO_MODE_ETSF if netcdf file, else IO_MODE_MPI if MPI supports it, finally IO_MODE_FORTRAN errmsg=String with error message. Use `if (len_trim(errmsg) /= 0) ABI_ERROR(errmsg)` to handle possible errors in the caller.
SOURCE
419 subroutine nctk_fort_or_ncfile(filename, iomode, errmsg) 420 421 character(len=*),intent(inout) :: filename 422 character(len=*),intent(out) :: errmsg 423 integer,intent(out) :: iomode 424 425 ! ********************************************************************* 426 errmsg = "" 427 428 ! Default value 429 #ifdef HAVE_MPI_IO 430 iomode = IO_MODE_MPI 431 #else 432 iomode = IO_MODE_FORTRAN 433 #endif 434 435 ! Checking the existence of data file 436 if (.not.file_exists(filename)) then 437 ! Trick needed to run Abinit test suite in netcdf mode. 438 if (file_exists(nctk_ncify(filename))) then 439 write(std_out,"(3a)")"- File: ",trim(filename)," does not exist but found netcdf file with similar name." 440 filename = nctk_ncify(filename); iomode = IO_MODE_ETSF 441 end if 442 if (.not. file_exists(filename)) then 443 errmsg = 'Missing file: '//trim(filename) 444 end if 445 end if 446 447 end subroutine nctk_fort_or_ncfile
m_nctk/nctk_get_dim [ Functions ]
[ Top ] [ m_nctk ] [ Functions ]
NAME
nctk_get_dim
FUNCTION
Get the value of a dimension from its name.
INPUTS
ncid=Netcdf identifier. dimname=Name of the dimension. [datamode]=If True, the nc file is set in data mode (default=False)
OUTPUT
dimlen=Value of the dimension.
SOURCE
1899 integer function nctk_get_dim(ncid, dimname, dimlen, datamode) result(ncerr) 1900 1901 !Arguments ------------------------------------ 1902 !scalars 1903 integer,intent(in) :: ncid 1904 character(len=*),intent(in) :: dimname 1905 integer,intent(out) :: dimlen 1906 logical,optional,intent(in) :: datamode 1907 1908 !Local variables------------------------------- 1909 !scalars 1910 integer :: dimid 1911 1912 ! ********************************************************************* 1913 1914 ncerr = nf90_noerr 1915 1916 if (present(datamode)) then 1917 if (datamode) then 1918 NCF_CHECK(nctk_set_datamode(ncid)) 1919 end if 1920 end if 1921 1922 NCF_CHECK(nf90_inq_dimid(ncid, dimname, dimid)) 1923 NCF_CHECK(nf90_inquire_dimension(ncid, dimid, len=dimlen)) 1924 1925 end function nctk_get_dim
m_nctk/nctk_idname [ Functions ]
[ Top ] [ m_nctk ] [ Functions ]
NAME
nctk_idname
FUNCTION
Return the nc identifier from the name of the variable
SOURCE
299 integer function nctk_idname(ncid, varname) result(varid) 300 301 !Arguments ------------------------------------ 302 integer,intent(in) :: ncid 303 character(len=*),intent(in) :: varname 304 305 !Local variables------------------------------- 306 !scalars 307 integer :: ncerr 308 character(len=1000) :: msg 309 310 ! ********************************************************************* 311 312 #ifdef HAVE_NETCDF 313 ncerr = nf90_inq_varid(ncid, varname, varid) 314 315 if (ncerr /= nf90_noerr) then 316 write(msg,'(6a)')& 317 "NetCDF library returned: `",trim(nf90_strerror(ncerr)), "`", ch10,& 318 "while trying to get the ncid of variable: ",trim(varname) 319 ABI_ERROR(msg) 320 end if 321 #else 322 ABI_ERROR("Netcdf support is not activated") 323 write(std_out,*)ncid,varname 324 #endif 325 326 end function nctk_idname
m_nctk/nctk_ncify [ Functions ]
[ Top ] [ m_nctk ] [ Functions ]
NAME
nctk_ncify
FUNCTION
Append ".nc" to ipath if ipath does not end with ".nc"
SOURCE
340 function nctk_ncify(ipath) result(opath) 341 342 character(len=fnlen),intent(in) :: ipath 343 character(len=fnlen) :: opath 344 345 ! ********************************************************************* 346 347 if (.not. endswith(ipath, ".nc")) then 348 opath = trim(ipath)//'.nc' 349 else 350 opath = ipath 351 end if 352 353 end function nctk_ncify
m_nctk/nctk_open_create [ Functions ]
[ Top ] [ m_nctk ] [ Functions ]
NAME
nctk_open_create
FUNCTION
Create and open the netcdf file. Return exis status.
INPUTS
path=Name of the file comm=MPI communicator.
OUTPUT
ncid=Netcdf identifier.
SOURCE
736 integer function nctk_open_create(ncid, path, comm) result(ncerr) 737 738 !Arguments ------------------------------------ 739 integer,intent(out) :: ncid 740 integer,intent(in) :: comm 741 character(len=*),intent(in) :: path 742 743 !Local variables------------------------------- 744 integer :: input_len, cmode 745 character(len=strlen) :: my_string 746 747 ! ********************************************************************* 748 749 ! Always use mpiio mode (i.e. hdf5) if available so that one can perform parallel IO 750 if (nctk_has_mpiio) then 751 ncerr = nf90_einval 752 #ifdef HAVE_NETCDF_MPI 753 call wrtout(std_out, sjoin("- Creating HDf5 file with MPI-IO support:", path)) 754 ! Believe it or not, I have to use xmpi_comm_self even in sequential to avoid weird SIGSEV in the MPI layer! 755 ncerr = nf90_create(path, cmode=ior(ior(nf90_netcdf4, nf90_mpiio), nf90_write), ncid=ncid, & 756 comm=comm, info=xmpio_info) 757 #endif 758 else 759 ! Note that here we don't enforce nf90_netcdf4 hence the netcdf file with be in classic model. 760 call wrtout(std_out, sjoin("- Creating netcdf file WITHOUT MPI-IO support:", path)) 761 !ncerr = nf90_create(path, ior(nf90_clobber, nf90_write), ncid) 762 cmode = def_cmode_for_seq_create 763 ncerr = nf90_create(path, cmode=cmode, ncid=ncid) 764 if (xmpi_comm_size(comm) > 1) then 765 ABI_WARNING("netcdf without MPI-IO support with nprocs > 1!") 766 end if 767 end if 768 NCF_CHECK(ncerr) 769 770 ! Write etsf_header: file format, version and conventions. 771 NCF_CHECK(nf90_put_att(ncid, NF90_GLOBAL, "file_format", etsfio_file_format)) 772 NCF_CHECK(nf90_put_att(ncid, NF90_GLOBAL, "file_format_version", etsfio_version)) 773 NCF_CHECK(nf90_put_att(ncid, NF90_GLOBAL, "Conventions", etsfio_conventions)) 774 775 ! Add info on the code that produced this file. These are extensions not in the standard. 776 NCF_CHECK(nf90_put_att(ncid, NF90_GLOBAL, "code", "Abinit")) 777 NCF_CHECK(nf90_put_att(ncid, NF90_GLOBAL, "abinit_version", abinit_version)) 778 779 ! Define the basic dimensions used in ETSF-IO files. 780 NCF_CHECK(nctk_def_basedims(ncid, defmode=.True.)) 781 782 if (len_trim(INPUT_STRING) /= 0) then 783 ! Write string with input. 784 my_string = trim(INPUT_STRING) 785 if (DTSET_IDX /= -1 .and. index(INPUT_STRING, "jdtset ") == 0) then 786 my_string = "jdtset " // trim(itoa(DTSET_IDX)) // " " // trim(INPUT_STRING) 787 end if 788 789 input_len = len_trim(my_string) 790 NCF_CHECK(nctk_def_dims(ncid, nctkdim_t("input_len", input_len))) 791 NCF_CHECK(nctk_def_arrays(ncid, nctkarr_t("input_string", "c", "input_len"))) 792 !print *, trim(INPUT_STRING) 793 794 if (xmpi_comm_rank(comm) == 0) then 795 NCF_CHECK(nctk_set_datamode(ncid)) 796 ! Pass my_string(1:input_len)) instead from trim(string) to avoid SIGSEV on higgs_intel_19.0_serial 797 NCF_CHECK(nf90_put_var(ncid, nctk_idname(ncid, "input_string"), my_string(1:input_len))) 798 NCF_CHECK(nctk_set_defmode(ncid)) 799 end if 800 end if 801 802 end function nctk_open_create
m_nctk/nctk_open_modify [ Functions ]
[ Top ] [ m_nctk ] [ Functions ]
NAME
nctk_open_modfy
FUNCTION
Open an already existent netcdf file for modifications. Return exit status.
INPUTS
path=File name. comm=MPI communicator.
OUTPUT
ncid=Netcdf identifier.
SOURCE
823 integer function nctk_open_modify(ncid, path, comm) result(ncerr) 824 825 !Arguments ------------------------------------ 826 integer,intent(out) :: ncid 827 integer,intent(in) :: comm 828 character(len=*),intent(in) :: path 829 830 ! ********************************************************************* 831 832 if (.not. nctk_has_mpiio .and. xmpi_comm_size(comm) > 1) then 833 ABI_ERROR("netcdf without MPI-IO support and nprocs > 1!") 834 end if 835 836 if (xmpi_comm_size(comm) > 1 .or. nctk_has_mpiio) then 837 call wrtout(std_out, sjoin("- Opening HDf5 file with MPI-IO support:", path)) 838 #ifdef HAVE_NETCDF_MPI 839 ncerr = nf90_open_par(path, cmode=ior(ior(nf90_netcdf4, nf90_mpiio), nf90_write), & 840 comm=comm, info=xmpio_info, ncid=ncid) 841 NCF_CHECK_MSG(ncerr, sjoin("nf90_open_par: ", path)) 842 #else 843 ABI_ERROR("nprocs > 1 but netcdf does not support MPI-IO") 844 #endif 845 else 846 call wrtout(std_out, sjoin("- Opening netcdf file without MPI-IO support:", path)) 847 ncerr = nf90_open(path, nf90_write, ncid) 848 NCF_CHECK_MSG(ncerr, sjoin("nf90_open: ", path)) 849 end if 850 851 ! Set file in define mode. 852 NCF_CHECK(nctk_set_defmode(ncid)) 853 854 end function nctk_open_modify
m_nctk/nctk_open_read [ Functions ]
[ Top ] [ m_nctk ] [ Functions ]
NAME
nctk_open_read
FUNCTION
Open a netcdf file in read-only mode. Return exit status.
INPUTS
ncid=Netcdf identifier. comm=MPI communicator.
SOURCE
679 integer function nctk_open_read(ncid, path, comm) result(ncerr) 680 681 !Arguments ------------------------------------ 682 integer,intent(out) :: ncid 683 integer,intent(in) :: comm 684 character(len=*),intent(in) :: path 685 686 !Local variables------------------------------- 687 integer :: nprocs 688 689 ! ********************************************************************* 690 nprocs = xmpi_comm_size(comm) 691 692 ! Enforce netcdf4 only if the communicator contains more than one processor. 693 if (nctk_has_mpiio .and. nprocs > 1) then 694 #ifdef HAVE_NETCDF_MPI 695 ncerr = nf90_open(path, mode=ior(ior(nf90_netcdf4, nf90_mpiio), nf90_nowrite),& 696 comm=comm, info=xmpio_info, ncid=ncid) 697 #else 698 ncerr = nf90_einval 699 ABI_WARNING("Netcdf without MPI support. Cannot open file, will abort in caller") 700 #endif 701 NCF_CHECK_MSG(ncerr, sjoin("opening file:", path)) 702 else 703 ncerr = nf90_open(path, mode=nf90_nowrite, ncid=ncid) 704 NCF_CHECK_MSG(ncerr, sjoin("opening file:", path)) 705 !if (ncerr /= NC_EHDFERR) then 706 ! ncerr = nf90_open(path, mode=ior(ior(nf90_netcdf4), nf90_nowrite),& 707 ! comm=comm, info=xmpio_info, ncid=ncid) 708 !end if 709 if (nprocs > 1) then 710 ncerr = nf90_einval 711 ABI_WARNING("netcdf without MPI-IO support with nprocs > 1! Will abort in the caller") 712 end if 713 endif 714 715 end function nctk_open_read
m_nctk/nctk_read_datar [ Functions ]
[ Top ] [ m_nctk ] [ Functions ]
NAME
nctk_read_datar
FUNCTION
Read an array in real space in netcdf format
INPUTS
path=Filename varname=Name of the variable to read. ngfft(18)=information about 3D FFT cplex=1 for real arrays (e.g. GS rhor), 2 for complex array. nfft=number of points in the real space FFT mesh treated by this MPI proc nspden=number of spin-density components comm_fft=MPI communicator (used only if MPI-FFT). fftn3_distrib(n3)=rank of the processors which own fft planes in 3rd dimension. ffti3_local(n3)=local index for 3d dimension datar(cplex*nfft,nspden)= array in real space.
OUTPUT
Only writing
SOURCE
2157 integer function nctk_read_datar(path,varname,ngfft,cplex,nfft,nspden,& 2158 comm_fft,fftn3_distrib,ffti3_local,datar) result(ncerr) 2159 2160 !Arguments ------------------------------------ 2161 !scalars 2162 integer,intent(in) :: cplex,nfft,nspden,comm_fft 2163 character(len=*),intent(in) :: path,varname 2164 !arrays 2165 integer,intent(in) :: ngfft(18),fftn3_distrib(ngfft(3)),ffti3_local(ngfft(3)) 2166 real(dp),intent(out) :: datar(cplex*nfft,nspden) 2167 2168 !Local variables------------------------------- 2169 !scalars 2170 integer,parameter :: master=0 2171 integer :: ncid,varid,i3,nproc_fft,me_fft,i3_glob,n1,n2,n3,ispden 2172 logical :: ionode 2173 !arrays 2174 real(dp),allocatable :: glob_datar(:,:) 2175 2176 ! ************************************************************************* 2177 2178 nproc_fft = xmpi_comm_size(comm_fft); me_fft = xmpi_comm_rank(comm_fft) 2179 n1 = ngfft(1); n2 = ngfft(2); n3 = ngfft(3) 2180 2181 ! TODO: Be careful here because we should always create with HDF5 if available 2182 ! to avoid problems if we have to reread with nproc_fft > 1 and MPI-IO 2183 ionode = .True. 2184 if (nproc_fft == 1) then 2185 ncerr = nf90_open(path, mode=nf90_nowrite, ncid=ncid) 2186 else 2187 if (nctk_has_mpiio) then 2188 !write(std_out,*)"open_par: ",trim(path) 2189 !ncerr = nf90_open_par(path, nf90_nowrite, 2190 ! Don't know why but the format is not autodected! 2191 ncerr = nf90_einval 2192 #ifdef HAVE_NETCDF_MPI 2193 ncerr = nf90_open(path, mode=ior(ior(nf90_netcdf4, nf90_mpiio), nf90_nowrite),& 2194 comm=comm_fft, info=xmpio_info, ncid=ncid) 2195 #endif 2196 else 2197 ! MPI-FFT without MPI-support. Only master does IO 2198 ionode = (me_fft == master); ncerr = nf90_noerr 2199 if (ionode) ncerr = nf90_open(path, nf90_nowrite, ncid) 2200 end if 2201 end if 2202 NCF_CHECK_MSG(ncerr, sjoin("opening file: ",path)) 2203 2204 NCF_CHECK(nf90_inq_varid(ncid, varname, varid)) 2205 !write(std_out,*)"about to read varname, ngfft, cplex, nfft, nspden:", trim(varname), ngfft(:3), cplex,nfft,nspden 2206 2207 ! netcdf array has shape [cplex, n1, n2, n3, nspden] 2208 if (nproc_fft == 1) then 2209 ! No MPI-FFT --> easy 2210 NCF_CHECK(nf90_get_var(ncid, varid, datar, start=[1,1,1,1], count=[cplex, n1, n2, n3, nspden])) 2211 NCF_CHECK(nf90_close(ncid)) 2212 2213 else 2214 ! Handle data distribution. 2215 ABI_CHECK(mod(ngfft(3), nproc_fft) == 0, "assuming mod(n3, nproc_fft) == 0") 2216 2217 i3_glob = -1 2218 do i3=1,ngfft(3) 2219 if (fftn3_distrib(i3) == me_fft) then 2220 i3_glob = i3 2221 exit 2222 end if 2223 end do 2224 ABI_CHECK(i3_glob > 0, "negative i3_glob") 2225 2226 if (nctk_has_mpiio) then 2227 ! Use parallel IO with collective calls. 2228 ncerr = nf90_einval 2229 NCF_CHECK(nctk_set_collective(ncid, varid)) 2230 2231 do ispden=1,nspden 2232 ncerr = nf90_get_var(ncid, varid, datar(:,ispden), start=[1,1,1,i3_glob,ispden], & 2233 count=[cplex,ngfft(1),ngfft(2),ngfft(3)/nproc_fft,1]) 2234 NCF_CHECK(ncerr) 2235 end do 2236 else 2237 ! MPI-FFT without MPI-IO. Master read and broadcast (requires more memory and communication) 2238 ABI_MALLOC(glob_datar, (cplex*product(ngfft(1:3)), nspden)) 2239 if (ionode) then 2240 NCF_CHECK(nf90_get_var(ncid, varid, glob_datar, start=[1,1,1,1,1], count=[cplex,n1,n2,n3,nspden])) 2241 end if 2242 2243 call distrib_datar(ngfft,cplex,nfft,nspden,glob_datar,master,comm_fft,fftn3_distrib,ffti3_local,datar) 2244 ABI_FREE(glob_datar) 2245 end if 2246 2247 if (ionode) then 2248 NCF_CHECK(nf90_close(ncid)) 2249 end if 2250 end if 2251 2252 end function nctk_read_datar
m_nctk/nctk_set_atomic_units [ Functions ]
[ Top ] [ m_nctk ] [ Functions ]
NAME
nctk_set_atomic_units
FUNCTION
Set the attributes "units" to "atomic units" and "scale_to_atomic_units" to one.
INPUTS
ncid=Netcdf identifier. varname=Name of the variable
SOURCE
1189 integer function nctk_set_atomic_units(ncid, varname) result(ncerr) 1190 1191 !Arguments ------------------------------------ 1192 integer,intent(in) :: ncid 1193 character(len=*),intent(in) :: varname 1194 1195 !Local variables------------------------------- 1196 !scalars 1197 integer :: varid 1198 1199 ! ********************************************************************* 1200 1201 ncerr = nf90_noerr 1202 1203 varid = nctk_idname(ncid, varname) 1204 NCF_CHECK(nf90_put_att(ncid, varid, "units", "atomic units")) 1205 NCF_CHECK(nf90_put_att(ncid, varid, "scale_to_atomic_units", one)) 1206 1207 end function nctk_set_atomic_units
m_nctk/nctk_set_collective [ Functions ]
[ Top ] [ m_nctk ] [ Functions ]
NAME
nctk_set_collective
FUNCTION
Use collective IO for the given netcdf variable. Return exit status.
INPUTS
ncid=Netcdf file identifier. varid=Netcdf variable identifier.
SOURCE
1032 integer function nctk_set_collective(ncid, varid) result(ncerr) 1033 1034 !Arguments ------------------------------------ 1035 integer,intent(in) :: ncid,varid 1036 1037 ! *********************************************************************true 1038 1039 ncerr = nf90_einval 1040 #ifdef HAVE_NETCDF_MPI 1041 ncerr = nf90_var_par_access(ncid, varid, nf90_collective) 1042 #else 1043 ABI_ERROR("nctk_set_collective should not be called if NETCDF does not support MPI-IO") 1044 ABI_UNUSED((/ncid, varid/)) 1045 #endif 1046 1047 end function nctk_set_collective
m_nctk/nctk_set_datamode [ Functions ]
[ Top ] [ m_nctk ] [ Functions ]
NAME
nctk_set_datamode
FUNCTION
Set the file in data mode. Return exit status
INPUTS
ncid=Netcdf identifier. [reserve]
OUTPUT
ncerr=Exit status
SOURCE
980 integer function nctk_set_datamode(ncid, reserve) result(ncerr) 981 982 !Arguments ------------------------------------ 983 integer,intent(in) :: ncid 984 logical,intent(in),optional :: reserve 985 986 !Local variables------------------------------- 987 !scalars 988 logical :: do_reserve 989 990 ! ********************************************************************* 991 992 do_reserve = .False.; if (present(reserve)) do_reserve = reserve 993 994 ncerr = nf90_enddef(ncid) 995 996 ! Use same trick as in etsf_io 997 ! neded otherwise netcdf complains if the file is already in def mode. 998 if (ncerr /= nf90_noerr .and. ncerr /= -38) then 999 NCF_CHECK(ncerr) 1000 else 1001 ncerr = nf90_noerr 1002 end if 1003 1004 return 1005 1006 ! TODO 1007 if (do_reserve) then 1008 ncerr = nf90_enddef(ncid) 1009 !ncerr = nf90__enddef(ncid) 1010 else 1011 ncerr = nf90_enddef(ncid) 1012 end if 1013 1014 end function nctk_set_datamode
m_nctk/nctk_set_default_for_seq [ Functions ]
[ Top ] [ m_nctk ] [ Functions ]
NAME
nctk_set_default_for_seq
FUNCTION
Use netcdf classic mode for new files when only sequential-IO needs to be performed
SOURCE
277 subroutine nctk_use_classic_for_seq() 278 279 ! ********************************************************************* 280 281 #ifdef HAVE_NETCDF 282 ! Use netcdf classic mode. 283 def_cmode_for_seq_create = ior(nf90_clobber, nf90_write) 284 ABI_COMMENT("Using netcdf-classic mode") 285 #endif 286 287 end subroutine nctk_use_classic_for_seq
m_nctk/nctk_set_defmode [ Functions ]
[ Top ] [ m_nctk ] [ Functions ]
NAME
nctk_set_defmode
FUNCTION
Set the file in define mode, return exit status.
INPUTS
ncid=Netcdf identifier.
SOURCE
944 integer function nctk_set_defmode(ncid) result(ncerr) 945 946 !Arguments ------------------------------------ 947 integer,intent(in) :: ncid 948 949 ! ********************************************************************* 950 951 ncerr = nf90_redef(ncid) 952 ! Use same trick as in etsf_io 953 if (ncerr /= nf90_noerr .and. ncerr /= -39) then 954 NCF_CHECK(ncerr) 955 else 956 ncerr = nf90_noerr 957 end if 958 959 end function nctk_set_defmode
m_nctk/nctk_string_from_occopt [ Functions ]
[ Top ] [ m_nctk ] [ Functions ]
NAME
nctk_string_from_occopt
FUNCTION
SOURCE
367 pure function nctk_string_from_occopt(occopt) result(smearing) 368 369 integer,intent(in) :: occopt 370 character(len=etsfio_charlen) :: smearing 371 372 ! ********************************************************************* 373 374 select case (occopt) 375 case (3) 376 smearing = "Fermi-Dirac" 377 case (4) 378 smearing = "cold smearing of N. Marzari with minimization of the bump" 379 case (5) 380 smearing = "cold smearing of N. Marzari with monotonic function in the tail" 381 case (6) 382 smearing = "Methfessel and Paxton" 383 case (7) 384 smearing = "gaussian" 385 case (8) 386 smearing = "uniform" 387 case default 388 smearing = "none" 389 end select 390 391 end function nctk_string_from_occopt
m_nctk/nctk_test_mpiio [ Functions ]
[ Top ] [ m_nctk ] [ Functions ]
NAME
nctk_test_mpiio
FUNCTION
Test at run-time whether the netcdf library supports parallel IO and set the value of the module variable `nctk_has_mpiio`. This is a COLLECTIVE routine that should be called by all processors in MPI_COMM_WORLD at the beginning of the calculation
INPUTS
[print_warning]=TRUE if a warning about paral_kgb use has to be printed Optional, default=yes
SOURCE
519 subroutine nctk_test_mpiio(print_warning) 520 521 logical,intent(in),optional :: print_warning 522 523 !Local variables------------------------------- 524 !scalars 525 logical :: my_print_warning 526 #ifdef HAVE_NETCDF_MPI 527 integer,parameter :: master=0 528 integer :: ierr,ncid,ncerr 529 character(len=500) :: msg 530 character(len=fnlen) :: apath 531 #endif 532 533 ! ********************************************************************* 534 535 nctk_has_mpiio = .False. 536 my_print_warning=.true. ; if (present(print_warning)) my_print_warning=print_warning 537 538 #ifdef HAVE_NETCDF_MPI 539 if (xmpi_comm_rank(xmpi_world) == master) then 540 ! Try to open a file with hdf5. 541 apath = pick_aname() 542 ncerr = nf90_create(apath, cmode=ior(ior(nf90_netcdf4, nf90_mpiio), nf90_write), ncid=ncid, & 543 comm=xmpi_comm_self, info=xmpio_info) 544 545 if (ncerr == nf90_noerr) then 546 nctk_has_mpiio = .True. 547 call wrtout(std_out," Netcdf library supports MPI-IO", "COLL") 548 else if (ncerr == nf90_enopar) then 549 ! This is the value returned by the C function ifndef USE_PARALLEL 550 ABI_WARNING(sjoin("Netcdf lib does not support MPI-IO and: ", nf90_strerror(ncerr))) 551 nctk_has_mpiio = .False. 552 else 553 ! Maybe something wrong in the low-level layer! 554 ABI_WARNING(sjoin("Strange, netcdf seems to support MPI-IO but: ", nf90_strerror(ncerr))) 555 nctk_has_mpiio = .False. 556 end if 557 558 ncerr = nf90_close(ncid) 559 call delete_file(apath, ierr) 560 end if 561 562 ! Master broadcast nctk_has_mpiio 563 call xmpi_bcast(nctk_has_mpiio,master,xmpi_world,ierr) 564 565 if ((.not.nctk_has_mpiio).and.my_print_warning) then 566 write(msg,"(5a)") & 567 "The netcdf library does not support parallel IO, see message above",ch10,& 568 "Abinit won't be able to produce files in parallel e.g. when paral_kgb==1 is used.",ch10,& 569 "Action: install a netcdf4+HDF5 library with MPI-IO support." 570 ABI_WARNING(msg) 571 end if 572 #endif 573 574 #ifdef HAVE_NETCDF_DEFAULT 575 if (.not. nctk_has_mpiio) then 576 ABI_ERROR("--netcdf-default is on but netcdf library does not support MPI-IO. Aborting now") 577 end if 578 #endif 579 580 end subroutine nctk_test_mpiio
m_nctk/nctk_try_fort_or_ncfile [ Functions ]
[ Top ] [ m_nctk ] [ Functions ]
NAME
nctk_try_fort_or_ncfile
FUNCTION
If filename does not exist, a similar file with extension `.nc` is tried This trick is used to run the Abinit test suite in netcdf mode without changing the input files. The modification (if any) is logged to unit (Default: std_out)
SIDE EFFECTS
filename=Tentative filename in input. Changed to netcdf file if input filename does not exist and a file with netcdf extension is found.
OUTPUT
errmsg=String with error message if return value /= 0
SOURCE
470 integer function nctk_try_fort_or_ncfile(filename, errmsg, unit) result(ierr) 471 472 character(len=*),intent(inout) :: filename 473 character(len=*),intent(out) :: errmsg 474 integer,optional,intent(in) :: unit 475 476 !Local variables------------------------------- 477 !scalars 478 integer :: unt 479 480 ! ********************************************************************* 481 482 unt = std_out; if (present(unit)) unt = unit 483 ierr = 0; errmsg = "" 484 485 if (.not.file_exists(filename)) then 486 ! Try netcdf exists. 487 if (file_exists(nctk_ncify(filename))) then 488 if (unt /= dev_null) then 489 write(unt,"(3a)")"- File: ",trim(filename)," does not exist but found netcdf file with similar name." 490 end if 491 filename = nctk_ncify(filename) 492 end if 493 if (.not. file_exists(filename)) then 494 ierr = 1; errmsg = 'Cannot find file: '//trim(filename) 495 end if 496 end if 497 498 end function nctk_try_fort_or_ncfile
m_nctk/nctk_write_datar [ Functions ]
[ Top ] [ m_nctk ] [ Functions ]
NAME
nctk_write_datar
FUNCTION
Write an array in real space in netcdf format
INPUTS
path=Filename varname=Name of the variable to write. ngfft(18)=information about 3D FFT cplex=1 for real arrays (e.g. GS rhor), 2 for complex array. nfft=number of points in the real space FFT mesh treated by this MPI proc nspden=number of spin-density components comm_fft=MPI communicator (used only if MPI-FFT). fftn3_distrib(n3)=rank of the processors which own fft planes in 3rd dimension. ffti3_local(n3)=local index for 3d dimension datar(cplex*nfft,nspden)= array in real space. [action]
OUTPUT
Only writing
SOURCE
1955 integer function nctk_write_datar(varname,path,ngfft,cplex,nfft,nspden,& 1956 comm_fft,fftn3_distrib,ffti3_local,datar,action) result(ncerr) 1957 1958 !Arguments ------------------------------------ 1959 !scalars 1960 integer,intent(in) :: cplex,nfft,nspden,comm_fft 1961 character(len=*),intent(in) :: path,varname 1962 character(len=*),optional,intent(in) :: action 1963 !arrays 1964 integer,intent(in) :: ngfft(18),fftn3_distrib(ngfft(3)),ffti3_local(ngfft(3)) 1965 real(dp),target,intent(in) :: datar(cplex*nfft,nspden) 1966 1967 !Local variables------------------------------- 1968 !scalars 1969 integer,parameter :: master=0 1970 integer :: ncid,varid,i3,nproc_fft,me_fft,i3_glob,n1,n2,n3,ispden,cmode 1971 logical :: ionode 1972 character(len=nctk_slen) :: cplex_name,my_action 1973 !character(len=500) :: msg 1974 !arrays 1975 real(dp),allocatable :: glob_datar(:,:) 1976 1977 ! ************************************************************************* 1978 1979 ! FIXME: Default should be open but this enters into conflict with the abi_estf stuff! 1980 ! if we are using MPI-IO since file is not open with HDF5. 1981 !my_action = "create"; if (present(action)) my_action = action 1982 my_action = "open"; if (present(action)) my_action = action 1983 1984 nproc_fft = xmpi_comm_size(comm_fft); me_fft = xmpi_comm_rank(comm_fft) 1985 n1 = ngfft(1); n2 = ngfft(2); n3 = ngfft(3) 1986 1987 ! TODO: Be careful here because we should always create with HDF5 if available 1988 ! to avoid problems if we have to reread with nproc_fft > 1 and MPI-IO 1989 ionode = .True.; ncerr = nf90_noerr 1990 if (nproc_fft == 1) then 1991 select case(my_action) 1992 case ("open") 1993 ncerr = nf90_open(path, mode=nf90_write, ncid=ncid) 1994 case ("create") 1995 ncerr = nctk_open_create(ncid, path, comm_fft) 1996 case default 1997 ABI_ERROR(sjoin("Wrong action: ", my_action)) 1998 end select 1999 2000 else 2001 if (nctk_has_mpiio) then 2002 call wrtout(std_out, strcat("nctk_write_datar: using MPI-IO to write ", varname, path), "COLL") 2003 2004 ncerr = nf90_einval 2005 #ifdef HAVE_NETCDF_MPI 2006 select case(my_action) 2007 case ("open") 2008 ncerr = nf90_open(path, mode=nf90_write, comm=comm_fft, info=xmpio_info, ncid=ncid) 2009 case ("create") 2010 ncerr = nf90_create(path, cmode=ior(ior(nf90_netcdf4, nf90_mpiio), nf90_write), & 2011 comm=comm_fft, info=xmpio_info, ncid=ncid) 2012 case default 2013 ABI_ERROR(strcat("Wrong action:", my_action)) 2014 end select 2015 #endif 2016 else 2017 ! MPI-FFT without MPI-support. Only master performs IO 2018 ionode = (me_fft == master) 2019 if (ionode) then 2020 select case(my_action) 2021 case ("open") 2022 ncerr = nf90_open(path, mode=nf90_write, ncid=ncid) 2023 case ("create") 2024 !ncerr = nf90_create(path, cmode=nf90_clobber, ncid=ncid) 2025 cmode = def_cmode_for_seq_create 2026 ncerr = nf90_create(path, cmode=cmode, ncid=ncid) 2027 case default 2028 ABI_ERROR(strcat("Wrong action:", my_action)) 2029 end select 2030 end if 2031 end if 2032 end if 2033 NCF_CHECK_MSG(ncerr, sjoin("opening file:", path)) 2034 2035 if (ionode) then 2036 ! Define dims and variables 2037 !write(std_out,*)"defing dims",trim(varname)," in file: ",path 2038 cplex_name = strcat("real_or_complex_", varname) 2039 ncerr = nctk_def_dims(ncid, [& 2040 nctkdim_t(cplex_name, cplex),& 2041 nctkdim_t("number_of_grid_points_vector1", n1),& 2042 nctkdim_t("number_of_grid_points_vector2", n2),& 2043 nctkdim_t("number_of_grid_points_vector3", n3),& 2044 nctkdim_t("number_of_components", nspden)], defmode=.True.) 2045 NCF_CHECK(ncerr) 2046 2047 ncerr = nctk_def_one_array(ncid, nctkarr_t(name=varname, dtype="dp", shape_str=strcat(cplex_name, & 2048 ", number_of_grid_points_vector1, number_of_grid_points_vector2, number_of_grid_points_vector3, number_of_components")),& 2049 varid=varid) 2050 2051 ! Add attributes 2052 varid = nctk_idname(ncid, varname) 2053 NCF_CHECK(nf90_put_att(ncid, varid, "units", "atomic units")) 2054 NCF_CHECK(nf90_put_att(ncid, varid, "scale_to_atomic_units", one)) 2055 end if 2056 2057 if (nproc_fft == 1) then 2058 ! no MPI-FFT --> write data directly. 2059 varid = nctk_idname(ncid, varname) 2060 NCF_CHECK(nctk_set_datamode(ncid)) 2061 NCF_CHECK(nf90_put_var(ncid, varid, datar, start=[1,1,1,1,1], count=[cplex, n1, n2, n3, nspden])) 2062 NCF_CHECK(nf90_close(ncid)) 2063 2064 else 2065 ! Must handle data distribution. 2066 ABI_CHECK(mod(n3, nproc_fft) == 0, "assuming mod(n3, nproc_fft) == 0") 2067 2068 i3_glob = -1 2069 do i3=1,ngfft(3) 2070 if (fftn3_distrib(i3) == me_fft) then 2071 i3_glob = i3 2072 exit 2073 end if 2074 end do 2075 ABI_CHECK(i3_glob > 0, "negative i3_glob") 2076 2077 !do i3=i3_glob,ngfft(3) 2078 ! if (fftn3_distrib(i3) /= me_fft) exit 2079 ! !assert all(ffn3_distrib(i3_glob:i3_glob -1 + ngfft(3) / nproc_fft) == me_fft) 2080 !end do 2081 !i3_glob = i3- 1 2082 !print*,"i3_glob",i3_glob 2083 2084 if (ionode) then 2085 NCF_CHECK(nf90_enddef(ncid)) 2086 end if 2087 2088 ! Array on disk has shape [cplex, n1, n2, n3, nspden] 2089 if (nctk_has_mpiio) then 2090 ! Use collective IO. 2091 ncerr = nf90_einval 2092 NCF_CHECK(nctk_set_collective(ncid, varid)) 2093 2094 do ispden=1,nspden 2095 ncerr = nf90_put_var(ncid, varid, datar(:,ispden), start=[1,1,1,i3_glob,ispden], & 2096 count=[cplex,n1,n2,n3/nproc_fft,1]) 2097 NCF_CHECK(ncerr) 2098 end do 2099 else 2100 ! MPI-FFT without MPI-IO. Collect data (requires more memory and communication) 2101 ABI_MALLOC(glob_datar, (cplex*product(ngfft(1:3)), nspden)) 2102 call collect_datar(ngfft,cplex,nfft,nspden,datar,comm_fft,fftn3_distrib,ffti3_local,glob_datar,master=master) 2103 2104 if (ionode) then 2105 ! Write global array. 2106 NCF_CHECK(nf90_put_var(ncid, varid, glob_datar, start=[1,1,1,1,1], count=[cplex,n1,n2,n3,nspden])) 2107 end if 2108 ABI_FREE(glob_datar) 2109 end if 2110 2111 if (ionode) then 2112 NCF_CHECK(nf90_close(ncid)) 2113 end if 2114 end if 2115 2116 !ok = .True. 2117 ! Sequential IO 2118 !do rank=0,nproc_fft-1 2119 ! if (rank == me_fft) then 2120 ! ncerr = nf90_open(path, mode=nf90_write, ncid=ncid) 2121 ! do ispden=1,nspden 2122 ! ncerr = nf90_put_var(ncid, varid, datar(1:,ispden), start=[1,1,1,i3_glob,ispden], & 2123 ! count=[cplex,ngfft(1),ngfft(2),ngfft(3)/nproc_fft,1]) 2124 ! end do 2125 ! ncerr = nf90_close(ncid) 2126 ! end if 2127 ! call xmpi_barrier(comm_fft) 2128 !end do 2129 2130 end function nctk_write_datar
m_nctk/nctk_write_dpscalars [ Functions ]
[ Top ] [ m_nctk ] [ Functions ]
NAME
nctk_write_dpscalars
FUNCTION
Write a list of **scalar** real(dp) variables. Return immediately if error
INPUTS
ncid=Netcdf identifier. varnames(:)=List of strings with the name of the variables values(:)=List of real(dp) scalars [datamode]=If True, the nc file is set in data mode (default=False)
SOURCE
1699 integer function nctk_write_dpscalars(ncid, varnames, values, datamode) result(ncerr) 1700 1701 !Arguments ------------------------------------ 1702 !scalars 1703 integer,intent(in) :: ncid 1704 logical,optional,intent(in) :: datamode 1705 !arrays 1706 real(dp),intent(in) :: values(:) 1707 character(len=*),intent(in) :: varnames(:) 1708 1709 !Local variables------------------------------- 1710 !scalars 1711 integer :: ii,varid 1712 1713 ! ********************************************************************* 1714 1715 ncerr = nf90_noerr 1716 1717 ABI_CHECK(size(varnames) == size(values), "Different size in varnames, values") 1718 1719 if (present(datamode)) then 1720 if (datamode) then 1721 NCF_CHECK(nctk_set_datamode(ncid)) 1722 end if 1723 end if 1724 1725 do ii=1,size(varnames) 1726 NCF_CHECK(nf90_inq_varid(ncid, varnames(ii), varid)) 1727 NCF_CHECK(nf90_put_var(ncid, varid, values(ii))) 1728 end do 1729 1730 end function nctk_write_dpscalars
m_nctk/nctk_write_ibz [ Functions ]
[ Top ] [ m_nctk ] [ Functions ]
NAME
nctk_write_ibz
FUNCTION
Write the list of the k-points in the IBZ with the corresponding weights in Netcdf format. Mainly used for passing data to AbiPy. This routine should be called by master only.
INPUTS
fname=File name kpoints(:,:)=List of k-points weights(:)=K-point weights
OUTPUT
ncerr=Exit status
SOURCE
1842 integer function nctk_write_ibz(fname, kpoints, weights) result(ncerr) 1843 1844 !Arguments ------------------------------------ 1845 !scalars 1846 character(len=*),intent(in) :: fname 1847 !arrays 1848 real(dp),intent(in) :: kpoints(:,:),weights(:) 1849 1850 !Local variables------------------------------- 1851 !scalars 1852 integer :: nkpts,ncid 1853 1854 ! ********************************************************************* 1855 1856 ABI_CHECK(size(kpoints, dim=2) == size(weights), "size(kpoints, dim=2) != size(weights)") 1857 nkpts = size(kpoints, dim=2) 1858 1859 NCF_CHECK_MSG(nctk_open_create(ncid, fname, xmpi_comm_self), sjoin("Creating:", fname)) 1860 1861 ncerr = nctk_def_dims(ncid, [ & 1862 nctkdim_t("number_of_reduced_dimensions",3), nctkdim_t("number_of_kpoints", nkpts)], defmode=.True.) 1863 NCF_CHECK(ncerr) 1864 1865 ncerr = nctk_def_array_list(ncid, [& 1866 nctkarr_t('reduced_coordinates_of_kpoints', "dp", "number_of_reduced_dimensions, number_of_kpoints"),& 1867 nctkarr_t('kpoint_weights', "dp", "number_of_kpoints")]) 1868 NCF_CHECK(ncerr) 1869 1870 NCF_CHECK(nctk_set_datamode(ncid)) 1871 1872 ncerr = nf90_put_var(ncid, nctk_idname(ncid, 'reduced_coordinates_of_kpoints'), kpoints) 1873 NCF_CHECK(ncerr) 1874 1875 NCF_CHECK(nf90_put_var(ncid, nctk_idname(ncid, 'kpoint_weights'), weights)) 1876 1877 NCF_CHECK(nf90_close(ncid)) 1878 1879 end function nctk_write_ibz
m_nctk/nctk_write_iscalars [ Functions ]
[ Top ] [ m_nctk ] [ Functions ]
NAME
nctk_write_iscalars
FUNCTION
Write a list of **scalar** integer variables. Return immediately if error
INPUTS
ncid=Netcdf identifier. varnames(:)=List of strings with the name of the variables values(:)=List of integer scalars [datamode]=If True, the nc file is set in data mode (default=False)
OUTPUT
ncerr=Exit status
SOURCE
1651 integer function nctk_write_iscalars(ncid, varnames, values, datamode) result(ncerr) 1652 1653 !Arguments ------------------------------------ 1654 !scalars 1655 integer,intent(in) :: ncid 1656 logical,optional,intent(in) :: datamode 1657 !arrays 1658 integer,intent(in) :: values(:) 1659 character(len=*),intent(in) :: varnames(:) 1660 1661 !Local variables------------------------------- 1662 !scalars 1663 integer :: ii,varid 1664 1665 ! ********************************************************************* 1666 1667 ABI_CHECK(size(varnames) == size(values), "Different size in varnames, values") 1668 1669 ncerr = nf90_noerr 1670 if (present(datamode)) then 1671 if (datamode) then 1672 NCF_CHECK(nctk_set_datamode(ncid)) 1673 end if 1674 end if 1675 1676 do ii=1,size(varnames) 1677 NCF_CHECK_MSG(nf90_inq_varid(ncid, varnames(ii), varid), sjoin("Inquiring: ", varnames(ii))) 1678 NCF_CHECK(nf90_put_var(ncid, varid, values(ii))) 1679 end do 1680 1681 end function nctk_write_iscalars
m_nctk/nctkarr_t [ Types ]
NAME
nctkarr_t
FUNCTION
Stores the name and the shape of a netcdf array
SOURCE
135 type,public :: nctkarr_t 136 character(len=nctk_slen) :: name ! name of the array. 137 character(len=4) :: dtype ! string specifying the type. 138 character(len=nctk_slen) :: shape_str ! string with the shape. e.g. "dim1, dim2" for [dim1, dim2] array. 139 end type nctkarr_t
m_nctk/nctkerr_t [ Types ]
NAME
nctkerr_t
FUNCTION
SOURCE
96 type,private :: nctkerr_t 97 #ifdef HAVE_NETCDF 98 integer :: ncerr = nf90_noerr 99 #else 100 integer :: ncerr = 0 101 #endif 102 integer :: line = 0 103 character(len=fnlen) :: file = "Dummy File" 104 character(len=2048) :: msg="No error detected" 105 end type nctkerr_t
m_nctk/nctkvar_t [ Structures ]
[ Top ] [ m_nctk ] [ Structures ]
NAME
nctkvar_t
FUNCTION
This structure stores variable information, such as name, NetCDF id, type, shape and dimensions. It contains the following elements:
SOURCE
152 type nctkvar_t 153 154 integer :: id 155 ! the id used by NetCDF to access this variable. 156 157 integer :: xtype 158 ! the type of the variable 159 160 integer :: ndims 161 ! the number of dimensions (0 for scalar variable). 162 163 integer :: natts 164 ! The number of attributes associated to the variable 165 166 character(len=nctk_slen) :: name 167 ! the variable name. 168 169 character(len=nctk_slen) :: dimnames(nctk_max_dims) 170 ! the name corresponding to each dimension 171 ! Only if array variable, only (1:ndims) are relevent 172 173 integer :: dimids(nctk_max_dims) = -1 174 ! The id of the dimensions. only (1:ndims) are relevent 175 176 integer :: dimlens(nctk_max_dims) = 0 177 ! the size for each dimension if array variable, only (1:ndims) are relevent 178 179 !character(len=nctk_slen) :: dimnames(nctk_max_dims) 180 !character(len=nctk_slen), pointer :: ncattnames(:) 181 ! * ncattnames: the name corresponding to all associated attributes 182 183 end type nctkvar_t
m_nctk/str2xtype [ Functions ]
[ Top ] [ m_nctk ] [ Functions ]
NAME
str2xtype
FUNCTION
Return the netcdf type from a string. Possible values: c or ch for NF90_CHAR i or int for NF90_INTtrue sp for NF90_FLOAT dp for NF90_DOUBLE
SOURCE
598 integer function str2xtype(string) result(xtype) 599 600 !Arguments ------------------------------------ 601 character(len=*),intent(in) :: string 602 603 ! ********************************************************************* 604 605 !Type FORTRAN API Mnemonic Bits 606 !byte NF90_BYTE 8 607 !char NF90_CHAR 8 608 !short NF90_SHORT 16 609 !int NF90_INT 32 610 !float NF90_FLOAT 32 611 !double NF90_DOUBLE 64 612 613 select case (string) 614 case ("c", "ch", "char") 615 xtype = NF90_CHAR 616 case ("i", "int") 617 xtype = NF90_INT 618 case ("sp") 619 xtype = NF90_FLOAT 620 case ("dp") 621 xtype = NF90_DOUBLE 622 case default 623 ABI_ERROR(sjoin("Invalid string type:", string)) 624 end select 625 626 end function str2xtype
m_nctk/var_from_id [ Functions ]
[ Top ] [ m_nctk ] [ Functions ]
NAME
var_from_id
FUNCTION
Initialize a nctkvar_t object from the variable id
INPUTS
ncid=NC file handle varid=Variable ID
OUTPUT
var<nctkvar_t>=Info on the variable.
SOURCE
2419 subroutine var_from_id(ncid, varid, var) 2420 2421 !Arguments ------------------------------------ 2422 integer, intent(in) :: ncid, varid 2423 type(nctkvar_t), intent(out) :: var 2424 2425 !Local variables------------------------------- 2426 !scalars 2427 integer :: ii, ncerr 2428 !character(len=NF90_MAX_NAME) :: ncname 2429 2430 ! ********************************************************************* 2431 2432 ! Get info about the variable. 2433 var%id = varid 2434 ncerr = nf90_inquire_variable(ncid, var%id, & 2435 name=var%name, xtype=var%xtype, ndims=var%ndims, dimids=var%dimids, natts=var%natts) 2436 NCF_CHECK(ncerr) 2437 2438 ! Get info about dimensions. 2439 if (var%ndims > 0) then 2440 do ii=1,var%ndims 2441 ncerr = nf90_inquire_dimension(ncid, var%dimids(ii), len=var%dimlens(ii), name=var%dimnames(ii)) 2442 NCF_CHECK(ncerr) 2443 end do 2444 end if 2445 2446 ! Get the number of attributes and their names. 2447 !if (var%natts > 0) then 2448 ! do ii=1,var%natts 2449 ! ncerr = nf90_inq_attname(ncid, var%id, ii, var%attnames(ii)) 2450 ! end do 2451 !end if 2452 2453 end subroutine var_from_id
m_nctk/var_from_name [ Functions ]
[ Top ] [ m_nctk ] [ Functions ]
NAME
var_from_name
FUNCTION
Initialize a nctkvar_t object from the variable name
INPUTS
ncid=NC file handle name=Variable name
OUTPUT
var<nctkvar_t>=Info on the variable.
SOURCE
2474 subroutine var_from_name(ncid, name, var) 2475 2476 !Arguments ------------------------------------ 2477 integer, intent(in) :: ncid 2478 character(len=*),intent(in) :: name 2479 type(nctkvar_t), intent(out) :: var 2480 2481 !Local variables------------------------------- 2482 !scalars 2483 integer :: varid 2484 2485 ! ********************************************************************* 2486 2487 varid = nctk_idname(ncid, name) 2488 call var_from_id(ncid, varid, var) 2489 2490 end subroutine var_from_name
m_nctk/write_var_netcdf [ Functions ]
[ Top ] [ m_nctk ] [ Functions ]
NAME
write_var_netcdf
FUNCTION
Write variable into a netcdf dataset TODO: Remove!
INPUTS
arr_int arr_real marr narr typevar varname
OUTPUT
(only writing)
SOURCE
2746 subroutine write_var_netcdf(arr_int,arr_real,marr,narr,ncid,typevar,varname) 2747 2748 !Arguments ------------------------------------ 2749 !scalars 2750 integer,intent(in) :: narr,marr,ncid 2751 character(len=*),intent(in) :: varname 2752 character(len=3),intent(in) :: typevar 2753 !arrays 2754 integer,intent(in) :: arr_int(marr) 2755 real(dp),intent(in) :: arr_real(marr) 2756 2757 !Local variables------------------------------- 2758 !scalars 2759 integer :: var_id,var_type,vardim_id,ncerr 2760 !character(len=500) :: msg 2761 2762 ! ************************************************************************* 2763 2764 !write(std_out,*)"about to write varname: ",trim(varname) 2765 2766 #if defined HAVE_NETCDF 2767 if (ncid>0) then 2768 ! ### Put the file in definition mode 2769 ncerr=nf90_redef(ncid) 2770 if (ncerr/=NF90_NOERR.and.ncerr/=NF90_EINDEFINE) then 2771 NCF_CHECK_MSG(ncerr,'nf90_redef') 2772 end if 2773 ! ### Define the dimensions 2774 if (narr==1)then 2775 ncerr=nf90_inq_dimid(ncid,'one',vardim_id) 2776 NCF_CHECK_MSG(ncerr,'nf90_inq_varid') 2777 else 2778 ncerr=nf90_def_dim(ncid,trim(varname),narr,vardim_id) 2779 NCF_CHECK_MSG(ncerr,'nf90_def_dim') 2780 end if 2781 ! ### Define the variables 2782 if (typevar=='INT') then 2783 var_type=NF90_INT 2784 else if (typevar=='DPR') then 2785 var_type=NF90_DOUBLE 2786 end if 2787 ncerr=nf90_def_var(ncid, trim(varname), var_type, vardim_id, var_id) 2788 NCF_CHECK_MSG(ncerr,'nf90_def_var') 2789 ! ### Put the file in data mode 2790 ncerr=nf90_enddef(ncid) 2791 if (ncerr/=NF90_NOERR.and.ncerr/=NF90_ENOTINDEFINE) then 2792 NCF_CHECK_MSG(ncerr,'nf90_enddef') 2793 end if 2794 ! ### Write variables into the dataset 2795 if (typevar=='INT') then 2796 ncerr=nf90_put_var(ncid,var_id,arr_int,start=(/1/),count=(/narr/)) 2797 else if (typevar=='DPR') then 2798 ncerr=nf90_put_var(ncid,var_id,arr_real,start=(/1/),count=(/narr/)) 2799 end if 2800 NCF_CHECK_MSG(ncerr,'nf90_put_var') 2801 end if 2802 #endif 2803 2804 end subroutine write_var_netcdf