Skip to content
This repository

HTTPS clone URL

Subversion checkout URL

You can clone with HTTPS or Subversion.

Download ZIP
Fetching contributors…

Octocat-spinner-32-eaf2f5

Cannot retrieve contributors at this time

file 164 lines (129 sloc) 4.51 kb
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163
# Copyright (C) 2010,2011,2012 The ESPResSo project
#
# This file is part of ESPResSo.
#
# ESPResSo is free software: you can redistribute it and/or modify
# it under the terms of the GNU General Public License as published by
# the Free Software Foundation, either version 3 of the License, or
# (at your option) any later version.
#
# ESPResSo is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
# GNU General Public License for more details.
#
# You should have received a copy of the GNU General Public License
# along with this program. If not, see <http://www.gnu.org/licenses/>.

source "tests_common.tcl"

require_feature "DIPOLES"
require_feature "FFTW"
require_feature "ROTATION"
require_feature "CONSTRAINTS"
require_max_nodes_per_side 2

set tcl_precision 15

# volume fraction
set rho 0.10
set dipole_lambda 3.0
set particle_radius 1.0
set bjerrum 1
set int_time_step 0.002

set external_H_z 0.5

# this are NOT relative errors, note!
set accuracy_p3m 1.0e-3
set accuracy_mdlc 1.0e-8

set pos_dip_data [open "mdlc_system.data" r]

set counter 0
while { [eof $pos_dip_data] == 0 } {
    set line [gets $pos_dip_data]
    incr counter
}
close $pos_dip_data
set n_particle [expr $counter-1]

set dipole_modulus [expr sqrt($dipole_lambda * pow(2*$particle_radius,3))]

# setting box paramters
set box_l [expr sqrt($n_particle * 3.141592654 /$rho)*$particle_radius]
set skin 0.5

# give Espresso some parameters
setmd time_step $int_time_step
setmd skin $skin
setmd box_l $box_l $box_l $box_l
setmd periodic 1 1 1
setmd max_num_cells 2500

# read positions and torques from file
set pos_dip_data [open "mdlc_system.data" r]
for {set i 0} {$i < $n_particle} {incr i} {
    set line [gets $pos_dip_data]
    set partid [lindex $line 0]
    set posx [lindex $line 1]
    set posy [lindex $line 2]
    set posz [lindex $line 3]
    set posx [expr $posx ]
    set posy [expr $posy ]
    set posz [expr $posz ]
    
    set dipx [lindex $line 4]
    set dipy [lindex $line 5]
    set dipz [lindex $line 6]
    part $partid pos $posx $posy $posz dip $dipx $dipy $dipz
}
close $pos_dip_data

thermostat off

puts "MDLC test case is running ..."

#puts "\n Tuning p3m parameters...."
#puts [inter magnetic $bjerrum p3m tunev2 mesh 16 accuracy $accuracy_p3m]
#puts [inter magnetic]

# the following parameters where obtained by running the above tuning
inter magnetic $bjerrum p3m 26.7723205659024 16 4 0.109394851183892 9.97941572136326e-06
inter magnetic epsilon metallic

# calculating the gap-size
set gap_size [expr $box_l - (2*$particle_radius)]

#puts "\n switch on mdlc..."
inter magnetic mdlc $accuracy_mdlc $gap_size

integrate 0

set failed 0

set rms_force 0
set rms_torque 0

if { [catch {

    # compare energy
    set energy_data [open "mdlc_expected_energy.data" r]
    set line [gets $energy_data]
    set energy_from_mdds [lindex $line 0]
    close $energy_data
    set magnetic_energy [analyze energy magnetic]
    set energy_err [expr ($magnetic_energy - $energy_from_mdds)]

    # field is not considered in the energy
    constraint ext_magn_field 0 0 0

    if {$energy_err > $accuracy_p3m} {
set failed 1
puts "error of energy is $energy_err -> too large"
    }

    constraint ext_magn_field 0 0 $external_H_z

    # compare forces and torques with results from "Wang code"
    set force_torque_data [open "mdlc_expected_force_torque.data" r]
    for {set i 0} {$i < $n_particle} {incr i} {
set line [gets $force_torque_data]

set diff_force [vecsub [part $i print force] [lrange $line 1 3]]
set error_force [veclen $diff_force]
set rms_force [expr $rms_force + pow($error_force, 2)]

set diff_torque [vecsub [part $i print torque] [lrange $line 4 6]]
set error_torque [veclen $diff_torque]
set rms_torque [expr $rms_torque + pow($error_torque, 2)]
    }
    close $force_torque_data

    set rms_force [expr sqrt($rms_force/$n_particle)]
    if {$rms_force > $accuracy_p3m} {
puts "rms force error $rms_force -> too large"
set failed 1
    }

    set rms_torque [expr sqrt($rms_torque/$n_particle)]
    if {$rms_torque > 2*$accuracy_p3m} {
puts "rms torque error $rms_torque -> too large"
puts "we let it through for now, was always too high"
# set failed 1
    }

    if {$failed==1} {
error_exit "\nerror of forces, torques or magnetic energy was too large"
    }
} res ] } {
    error_exit $res
}

exit 0
Something went wrong with that request. Please try again.