TABLE OF CONTENTS
ABINIT/m_multipoles [ Modules ]
NAME
m_multipoles
FUNCTION
Compute spatial multipole moments of input array on FFT grid
COPYRIGHT
Copyright (C) 2003-2024 ABINIT group (MJV, MT, XG) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
SOURCE
16 #if defined HAVE_CONFIG_H 17 #include "config.h" 18 #endif 19 20 #include "abi_common.h" 21 22 module m_multipoles 23 24 use defs_basis 25 use m_errors 26 use m_abicore 27 use m_distribfft 28 use m_xmpi 29 use m_atomdata 30 31 use defs_abitypes, only : mpi_type 32 use m_io_tools, only : open_file 33 use m_geometry, only : xred2xcart 34 use m_mpinfo, only : ptabs_fourdp 35 36 implicit none 37 38 private
m_multipoles/multipoles_fftr [ Functions ]
[ Top ] [ m_multipoles ] [ Functions ]
NAME
multipoles_fftr
FUNCTION
Compute spatial multipole moments of input array on FFT grid Namely, the electrical dipole, quadrupole, etc... of the electron density
INPUTS
arraysp(nfft,nspden)=the array whose average has to be computed [distribfft<type(distribfft_type)>]= -optional- infos related to FFT parallelism [mpi_comm_grid]= -optional- MPI communicator over the grid origin(3)=vector defining the origin of the dipole (point of observation, in reduced coordinates) nfft=number of FFT points stored by one proc nfftot=total number of FFT points ngfft =number of subdivisions along each lattice vector nspden=number of spin-density components rprimd = dimensionful lattice vectors
OUTPUT
dipole(nspden)=mean value of the dipole of input array, for each nspden component
SOURCE
72 subroutine multipoles_fftr(arraysp,dipole,nfft,ngfft,nspden,rprimd,origin,& 73 distribfft,mpi_comm_grid) 74 75 !Arguments ------------------------------------ 76 !scalars 77 integer,intent(in) :: nfft,nspden 78 integer,intent(in),optional :: mpi_comm_grid 79 type(distribfft_type),intent(in),optional,target :: distribfft 80 !arrays 81 integer,intent(in) :: ngfft(18) 82 real(dp),intent(in) :: arraysp(nfft,nspden),origin(3),rprimd(3,3) 83 real(dp),intent(out) :: dipole(3,nspden) 84 !Local variables------------------------------- 85 !scalars 86 integer,parameter :: ishift=5 87 integer :: ierr,ifft1,ifft2,ifft3,ifft0,ifft,ispden,i1,i2,i3,n1,n2,n3 88 integer :: me_fft,my_mpi_comm,nfftot 89 logical :: fftgrid_found 90 real(dp) :: invn1,invn2,invn3,invnfftot,wrapfft1,wrapfft2,wrapfft3 91 character(len=500) :: msg 92 type(distribfft_type),pointer :: my_distribfft 93 !arrays 94 integer, ABI_CONTIGUOUS pointer :: fftn3_distrib(:),ffti3_local(:) 95 real(dp) :: dipole_tmp(3,nspden) 96 97 ! ************************************************************************* 98 99 100 !Several initializations 101 my_mpi_comm=xmpi_comm_self;if (present(mpi_comm_grid)) my_mpi_comm=mpi_comm_grid 102 n1=ngfft(1);n2=ngfft(2);n3=ngfft(3) 103 nfftot=n1*n2*n3;invnfftot=one/(dble(nfftot)) 104 invn1=one/real(n1,kind=dp);invn2=one/real(n2,kind=dp);invn3=one/real(n3,kind=dp) 105 me_fft=xmpi_comm_rank(my_mpi_comm) 106 107 dipole(:,:)=zero 108 109 !Get the distrib associated with the FFT grid 110 if (present(distribfft)) then 111 my_distribfft => distribfft 112 else 113 ABI_MALLOC(my_distribfft,) 114 call init_distribfft_seq(my_distribfft,'f',n2,n3,'fourdp') 115 end if 116 fftgrid_found=.false. 117 if (n2 == my_distribfft%n2_coarse ) then 118 if (n3 == size(my_distribfft%tab_fftdp3_distrib)) then 119 fftn3_distrib => my_distribfft%tab_fftdp3_distrib 120 ffti3_local => my_distribfft%tab_fftdp3_local 121 fftgrid_found=.true. 122 end if 123 end if 124 if (n2 == my_distribfft%n2_fine ) then 125 if (n3 == size(my_distribfft%tab_fftdp3dg_distrib)) then 126 fftn3_distrib => my_distribfft%tab_fftdp3dg_distrib 127 ffti3_local => my_distribfft%tab_fftdp3dg_local 128 fftgrid_found = .true. 129 end if 130 end if 131 if (.not.(fftgrid_found)) then 132 msg='Unable to find an allocated distrib for the FFT grid!' 133 ABI_BUG(msg) 134 end if 135 136 !Loop over FFT grid points 137 !$OMP PARALLEL PRIVATE(ifft,i1,i2,ifft1,ifft2,ispden,wrapfft1,wrapfft2) 138 do ifft3=1,n3 139 i3=mod(ifft3-1+ishift*n3,n3) 140 if(fftn3_distrib(1+i3)==me_fft) then 141 wrapfft3=mod(real(ifft3-1,kind=dp)*invn3-origin(3)+1.5_dp,one)-half 142 ifft0=1+n1*n2*(ffti3_local(1+i3)-1) 143 !$OMP SINGLE 144 dipole_tmp=zero 145 !$OMP END SINGLE 146 !$OMP DO COLLAPSE(2) REDUCTION(+:dipole_tmp) 147 do ifft2=1,n2 148 do ifft1=1,n1 149 i2=mod(ifft2-1+ishift*n2,n2) 150 i1=mod(ifft1-1+ishift*n1,n1) 151 wrapfft2=mod(real(ifft2-1,kind=dp)*invn2-origin(2)+1.5_dp,one)-half 152 wrapfft1=mod(real(ifft1-1,kind=dp)*invn1-origin(1)+1.5_dp,one)-half 153 ifft=ifft0+n1*i2+i1 154 155 ! Computation of integral(s) 156 do ispden=1,nspden 157 dipole_tmp(1,ispden)=dipole_tmp(1,ispden)+wrapfft1*arraysp(ifft,ispden) 158 dipole_tmp(2,ispden)=dipole_tmp(2,ispden)+wrapfft2*arraysp(ifft,ispden) 159 dipole_tmp(3,ispden)=dipole_tmp(3,ispden)+wrapfft3*arraysp(ifft,ispden) 160 end do 161 162 end do 163 end do 164 !$OMP END DO 165 !$OMP SINGLE 166 dipole=dipole+dipole_tmp 167 !$OMP END SINGLE 168 end if 169 end do 170 !$OMP END PARALLEL 171 172 !MPI parallelization 173 if (xmpi_comm_size(my_mpi_comm)>1) then 174 call xmpi_sum(dipole,my_mpi_comm,ierr) 175 end if 176 177 !From reduced to cartesian coordinates 178 do ispden=1,nspden 179 dipole(:,ispden)=matmul(rprimd,dipole(:,ispden))*invnfftot 180 end do 181 182 if (.not.present(distribfft)) then 183 call destroy_distribfft(my_distribfft) 184 ABI_FREE(my_distribfft) 185 end if 186 187 end subroutine multipoles_fftr
m_multipoles/multipoles_out [ Functions ]
[ Top ] [ m_multipoles ] [ Functions ]
NAME
multipoles_out
FUNCTION
Output multipole moments of input array on FFT grid, calculated with multipoles_fftr Namely, the electrical dipole, quadrupole, etc... of the electron density
INPUTS
rhor(nfft,nspden)=electronic density mpi_enreg=information about MPI parallelization natom=number of atoms nfft=number of FFT points stored by one proc ngfft =number of subdivisions along each lattice vector nspden=number of spin-density components ntypat=number of atom types rprimd(3,3)=dimensionful lattice vectors typat(ntypat)=types of atoms ucvol=unit cell volume unit_out=file unit to print out ziontypat(ntypat)=ionic charge of each type of atom
OUTPUT
NOTES
SOURCE
218 subroutine multipoles_out(rhor,mpi_enreg,natom,nfft,ngfft,nspden,& 219 ntypat,rprimd,typat,ucvol,unit_out,xred,ziontypat) 220 221 !Arguments ------------------------------------ 222 !scalars 223 integer,intent(in) :: natom,nfft,nspden,ntypat,unit_out 224 real(dp), intent(in) :: ucvol 225 type(MPI_type),intent(in) :: mpi_enreg 226 !arrays 227 integer, intent(in) :: ngfft(18),typat(natom) 228 real(dp),intent(in) :: rhor(nfft,nspden),rprimd(3,3),xred(3,natom),ziontypat(ntypat) 229 !Local variables ------------------------------ 230 !scalars 231 integer :: iatom,nspden_updn 232 real(dp) :: ziontotal 233 character(len=500) :: message 234 !arrays 235 real(dp) :: center_of_charge(3),dipole_el(3,2),dipole_ions_cart(3) 236 real(dp) :: dipole_ions_red(3),dipole_tot(3),tmp(3) 237 238 ! ************************************************************************* 239 240 !Separate spins only for nspden=2 241 nspden_updn=merge(1,2,nspden/=2) 242 243 !Title 244 245 !Get nuclear part of dipole 246 dipole_ions_red(:) = zero ; ziontotal = zero 247 do iatom = 1, natom 248 dipole_ions_red(:) = dipole_ions_red(:) + xred(:,iatom)*ziontypat(typat(iatom)) 249 ziontotal = ziontotal + ziontypat(typat(iatom)) 250 end do 251 dipole_ions_cart(:) = matmul(rprimd,dipole_ions_red(:)) 252 253 !Find coordinates of center of charge on FFT grid 254 center_of_charge(1:3) = dipole_ions_red(1:3)/ziontotal 255 256 !Get electronic part of dipole with respect to center of charge of ions (in cart. coord.) 257 dipole_el = zero 258 call multipoles_fftr(rhor(:,1:nspden_updn),dipole_el(:,1:nspden_updn),nfft,ngfft,nspden_updn,& 259 & rprimd,center_of_charge,distribfft=mpi_enreg%distribfft,mpi_comm_grid=mpi_enreg%comm_fft) 260 dipole_el(1:3,1:nspden_updn)=dipole_el(1:3,1:nspden_updn)*ucvol 261 262 !Take into account storage of rhor (up+dn,up) 263 if (nspden==2) then 264 tmp(1:3)=dipole_el(1:3,1) 265 dipole_el(1:3,1)=dipole_el(1:3,2) 266 dipole_el(1:3,2)=tmp(1:3)-dipole_el(1:3,2) 267 end if 268 269 !Compute total dipole 270 !NOTE: wrt center of charge, dipole_ions is 0 271 dipole_tot(1) = - sum(dipole_el(1,1:nspden_updn)) 272 dipole_tot(2) = - sum(dipole_el(2,1:nspden_updn)) 273 dipole_tot(3) = - sum(dipole_el(3,1:nspden_updn)) 274 275 !Output 276 write (message, '(2a)') ch10,' ----- Electric nuclear dipole wrt the center of ionic charge ----- ' 277 call wrtout(unit_out, message, 'COLL') 278 write (message, '(a,3(1x,ES12.5))') & 279 & ' Center of charge for ionic distribution (red. coord.): ',center_of_charge(1:3) 280 call wrtout(unit_out, message, 'COLL') 281 write (message, '(3a,3(1x,E16.6),3a,3(1x,E16.6),a)') ' -----',ch10,& 282 & ' Ionic dipole (cart. coord.) = ',dipole_ions_cart, ' (a.u.)',ch10, & 283 & ' = ',dipole_ions_cart/dipole_moment_debye,' (D)' 284 call wrtout(unit_out, message, 'COLL') 285 if (nspden/=2) then 286 !This is compatible with nspden=4 287 write (message, '(3a,3(1x,E16.6),3a,3(1x,E16.6),a)') ' -----',ch10,& 288 & ' Electronic dipole (cart. coord.)= ',dipole_el(:,1),' (a.u.)',ch10,& 289 & ' = ',dipole_el(:,1)/dipole_moment_debye,' (D)' 290 else 291 write (message, '(3a,3(1x,E16.6),a,3(2a,3(1x,E16.6),a))') ' -----',ch10,& 292 & ' Electronic dipole (cart. coord.)= ',dipole_el(:,1),' up (a.u.)',ch10,& 293 & ' = ',dipole_el(:,2),' dn (a.u.)',ch10,& 294 & ' = ',dipole_el(:,1)/dipole_moment_debye,' up (D)',ch10,& 295 & ' = ',dipole_el(:,2)/dipole_moment_debye,' dn (D)' 296 end if 297 call wrtout(unit_out, message, 'COLL') 298 write (message, '(3a,3(1x,E16.6),3a,3(1x,E16.6),a)') ' -----',ch10,& 299 & ' Total dipole (cart. coord.) = ',dipole_tot,' (a.u.)',ch10,& 300 & ' = ',dipole_tot/dipole_moment_debye,' (D)' 301 call wrtout(unit_out, message, 'COLL') 302 call wrtout(unit_out, ' ', 'COLL') 303 304 end subroutine multipoles_out
m_multipoles/out1dm [ Functions ]
[ Top ] [ m_multipoles ] [ Functions ]
NAME
out1dm
FUNCTION
Output the 1 dimensional mean of potential and density on the three reduced axis.
INPUTS
fnameabo_app_1dm=name of the file in which the data is written, appended with _1DM natom=number of atoms in unit cell mpi_enreg=information about MPI parallelization nfft=(effective) number of FFT grid points (for this processor) ngfft(18)=contain all needed information about 3D FFT, see ~abinit/doc/variables/vargs.htm#ngfft nspden=number of spin-density components ntypat=number of types of atoms in unit cell. rhor(nfft,nspden)=array for electron density in electrons/bohr**3. rprimd(3,3)=real space dimensional primitive translations (bohr) typat(natom)=type integer for each atom in cell ucvol=unit cell volume in bohr**3. vtrial(nfft,nspden)=INPUT Vtrial(r). xred(3,natom)=reduced coordinates of atoms znucl(ntypat)=real(dp), nuclear number of atom type
OUTPUT
SIDE EFFECTS
data written in file fnameabo_app_1dm
NOTES
SOURCE
339 subroutine out1dm(fnameabo_app_1dm,mpi_enreg,natom,nfft,ngfft,nspden,ntypat,& 340 rhor,rprimd,typat,ucvol,vtrial,xred,znucl) 341 342 !Arguments ------------------------------------ 343 !scalars 344 integer,intent(in) :: natom,nfft,nspden,ntypat 345 real(dp),intent(in) :: ucvol 346 character(len=fnlen),intent(in) :: fnameabo_app_1dm 347 type(MPI_type),intent(in) :: mpi_enreg 348 !arrays 349 integer,intent(in) :: ngfft(18),typat(natom) 350 real(dp),intent(in) :: rhor(nfft,nspden),rprimd(3,3),vtrial(nfft,nspden) 351 real(dp),intent(in) :: znucl(ntypat),xred(3,natom) 352 353 !Local variables------------------------------- 354 !scalars 355 integer :: ia,ib,idim,ierr,ifft,islice,ispden,na,nb,ndig,nslice,nu,temp_unit 356 integer :: comm_fft, nproc_fft, me_fft 357 real(dp) :: global_den,global_pot 358 character(len=2) :: symbol 359 character(len=500) :: message 360 type(atomdata_t) :: atom 361 !arrays 362 integer, ABI_CONTIGUOUS pointer :: fftn2_distrib(:),ffti2_local(:) 363 integer, ABI_CONTIGUOUS pointer :: fftn3_distrib(:),ffti3_local(:) 364 real(dp),allocatable :: lin_den(:),mean_pot(:),reduced_coord(:),xcart(:,:) 365 character(len=8),allocatable :: iden(:) 366 367 ! ************************************************************************* 368 369 !Initialize the file 370 write(message, '(a,a)' ) ' io1dm : about to open file ',trim(fnameabo_app_1dm) 371 call wrtout(std_out,message,'COLL') 372 call wrtout(ab_out,message,'COLL') 373 374 comm_fft = mpi_enreg%comm_fft; nproc_fft = xmpi_comm_size(comm_fft); me_fft = mpi_enreg%me_fft 375 376 if (me_fft == 0) then 377 if (open_file(fnameabo_app_1dm,message,newunit=temp_unit,status='unknown',form='formatted') /= 0) then 378 ABI_ERROR(message) 379 end if 380 rewind(temp_unit) 381 end if 382 383 write(message, '(a,a)' ) ch10,'# ABINIT package : 1DM file ' 384 if (me_fft == 0) write(temp_unit,'(a)') message 385 386 write(message, '(a,a)' )ch10,'# Primitive vectors of the periodic cell (bohr)' 387 if (me_fft == 0) write(temp_unit,'(a)') message 388 do nu=1,3 389 write(message, '(1x,a,i1,a,3f10.5)' ) '# R(',nu,')=',rprimd(:,nu) 390 if (me_fft == 0) write(temp_unit,'(a)') message 391 end do 392 393 write(message, '(a,a)' ) ch10,& 394 & '# Atom list Reduced coordinates Cartesian coordinates (bohr)' 395 if (me_fft == 0) write(temp_unit,'(a)') message 396 397 !Set up a list of character identifiers for all atoms : iden(ia) 398 ABI_MALLOC(iden,(natom)) 399 do ia=1,natom 400 call atomdata_from_znucl(atom,znucl(typat(ia))) 401 symbol = atom%symbol 402 403 ndig=int(log10(dble(ia)+0.5_dp))+1 404 if(ndig==1) write(iden(ia), '(a,a,i1,a)' )symbol,'(',ia,') ' 405 if(ndig==2) write(iden(ia), '(a,a,i1,a)' )symbol,'(',ia,') ' 406 if(ndig==3) write(iden(ia), '(a,a,i1,a)' )symbol,'(',ia,') ' 407 if(ndig==4) write(iden(ia), '(a,a,i1,a)' )symbol,'(',ia,')' 408 if(ndig>4)then 409 close(temp_unit) 410 write(message, '(a,i0)' )& 411 & ' out1dm : cannot handle more than 9999 atoms, while natom=',natom 412 ABI_ERROR(message) 413 end if 414 end do 415 416 !Compute cartesian coordinates, and print reduced and cartesian coordinates 417 ABI_MALLOC(xcart,(3,natom)) 418 call xred2xcart(natom,rprimd,xcart,xred) 419 do ia=1,natom 420 write(message, '(a,a,3f10.5,a,3f10.5)' ) & 421 & '# ',iden(ia),xred(1:3,ia),' ',xcart(1:3,ia) 422 if (me_fft == 0) write(temp_unit,'(a)') message 423 end do 424 ABI_FREE(iden) 425 ABI_FREE(xcart) 426 427 !Get the distrib associated with this fft_grid 428 call ptabs_fourdp(mpi_enreg,ngfft(2),ngfft(3),fftn2_distrib,ffti2_local,fftn3_distrib,ffti3_local) 429 430 do idim=1,3 431 432 nslice=ngfft(idim) 433 434 ! Dummy initialisation of na, nb and ifft. 435 na=0 ; nb=0; ifft=0 436 select case(idim) 437 case(1) 438 na=ngfft(2) ; nb=ngfft(3) 439 case(2) 440 na=ngfft(1) ; nb=ngfft(3) 441 case(3) 442 na=ngfft(1) ; nb=ngfft(2) 443 end select 444 445 ABI_MALLOC( reduced_coord,(nslice)) 446 ABI_MALLOC(mean_pot,(nslice)) 447 ABI_MALLOC(lin_den,(nslice)) 448 449 do ispden=1,nspden 450 451 if(ispden==1)then 452 write(message, '(a,a,a)' ) ch10,'#===========',& 453 & '=====================================================================' 454 if (me_fft == 0) write(temp_unit,'(a)') message 455 end if 456 457 select case(idim) 458 case(1) 459 write(message, '(a)' )'# Projection along the first dimension ' 460 case(2) 461 write(message, '(a)' )'# Projection along the second dimension ' 462 case(3) 463 write(message, '(a)' )'# Projection along the third dimension ' 464 end select 465 if (me_fft == 0) write(temp_unit,'(a)') message 466 467 if(nspden==2)then 468 select case(ispden) 469 case(1) 470 write(message, '(a)' )'# Spin up ' 471 case(2) 472 write(message, '(a)' )'# Spin down ' 473 end select 474 if (me_fft == 0) write(temp_unit,'(a)') message 475 else if (nspden == 4) then 476 select case(ispden) 477 case(1) 478 write(message, '(a)' )'# Spinor up up ' 479 case(2) 480 write(message, '(a)' )'# Spinor down down' 481 case(3) 482 write(message, '(a)' )'# Spinor up down' 483 case(4) 484 write(message, '(a)' )'# Spinor down up' 485 end select 486 if (me_fft == 0) write(temp_unit,'(a)') message 487 end if 488 489 write(message, '(2a)' ) ch10,& 490 & '# Red. coord. Mean KS potential Linear density ' 491 if (me_fft == 0) write(temp_unit,'(a)') message 492 493 write(message, '(a)' ) & 494 & '# (Hartree unit) (electron/red. unit)' 495 if (me_fft == 0) write(temp_unit,'(a)') message 496 497 global_pot=zero 498 global_den=zero 499 mean_pot(:)=zero 500 lin_den(:)=zero 501 do islice=1,nslice 502 reduced_coord(islice)=(islice-1)/dble(nslice) 503 end do 504 if (idim == 1) then 505 do islice=1,nslice 506 do ib=1,nb 507 ! skip z values not on this processor 508 if (fftn3_distrib(ib)/=mpi_enreg%me_fft) cycle 509 do ia=1,na 510 ifft = islice + nslice*( ia -1 + na *(ffti3_local(ib) -1) ) 511 mean_pot(islice)=mean_pot(islice)+vtrial(ifft,ispden) 512 lin_den(islice)=lin_den(islice)+rhor(ifft,ispden) 513 end do 514 end do 515 end do 516 else if (idim == 2) then 517 do islice=1,nslice 518 do ib=1,nb 519 ! skip z values not on this processor 520 if (fftn3_distrib(ib)/=mpi_enreg%me_fft) cycle 521 do ia=1,na 522 ifft = ia + na *( islice-1 + nslice*(ffti3_local(ib) -1) ) 523 mean_pot(islice)=mean_pot(islice)+vtrial(ifft,ispden) 524 lin_den(islice)=lin_den(islice)+rhor(ifft,ispden) 525 end do 526 end do 527 end do 528 else if (idim == 3) then 529 ! this full z slice is mine in parallel case 530 do islice=1,nslice 531 if (fftn3_distrib(islice)/=mpi_enreg%me_fft) cycle 532 do ib=1,nb 533 do ia=1,na 534 ifft = ia + na *( ib -1 + nb *(ffti3_local(islice)-1) ) 535 mean_pot(islice)=mean_pot(islice)+vtrial(ifft,ispden) 536 lin_den(islice)=lin_den(islice)+rhor(ifft,ispden) 537 end do 538 end do 539 end do 540 end if 541 mean_pot(:)=mean_pot(:)/dble(na*nb) 542 lin_den(:)=lin_den(:)/dble(na*nb)*ucvol 543 global_pot=sum(mean_pot)/dble(nslice) 544 global_den=sum(lin_den)/dble(nslice) 545 546 ! FFT parallelization 547 if(mpi_enreg%nproc_fft>1)then 548 call xmpi_sum(mean_pot ,mpi_enreg%comm_fft,ierr) 549 call xmpi_sum(global_pot,mpi_enreg%comm_fft,ierr) 550 call xmpi_sum(lin_den ,mpi_enreg%comm_fft,ierr) 551 call xmpi_sum(global_den,mpi_enreg%comm_fft,ierr) 552 end if 553 554 do islice=1,ngfft(idim) 555 write(message, '(f10.4,es20.6,es16.6)' )& 556 & reduced_coord(islice),mean_pot(islice),lin_den(islice) 557 if (me_fft == 0) write(temp_unit,'(a)') message 558 end do 559 560 write(message, '(a,a,es15.6,es16.6,a)' ) ch10,& 561 & '# Cell mean :',global_pot,global_den, ch10 562 if (me_fft == 0) write(temp_unit,'(a)') message 563 564 565 ! End of the loop on spins 566 end do 567 568 ABI_FREE(reduced_coord) 569 ABI_FREE(mean_pot) 570 ABI_FREE(lin_den) 571 572 ! End of the loops on the three dimensions 573 end do 574 575 if (me_fft == 0) then 576 write(message, '(a,a,a)' ) ch10,'#===========',& 577 & '=====================================================================' 578 call wrtout(temp_unit,message,'COLL') 579 close(temp_unit) 580 end if 581 582 end subroutine out1dm