~ubuntu-branches/debian/sid/lammps/sid

« back to all changes in this revision

Viewing changes to examples/USER/atc/elastic/in.electrostatic_bending

  • Committer: Package Import Robot
  • Author(s): Anton Gladky
  • Date: 2013-11-20 22:41:36 UTC
  • mfrom: (1.2.2)
  • Revision ID: package-import@ubuntu.com-20131120224136-tzx7leh606fqnckm
Tags: 0~20131119.git7162cf0-1
* [e65b919] Imported Upstream version 0~20131119.git7162cf0
* [f7bddd4] Fix some problems, introduced by upstream recently.
* [3616dfc] Use wrap-and-sort script.
* [7e92030] Ignore quilt dir

Show diffs side-by-side

added added

removed removed

Lines of Context:
 
1
echo both
 
2
units           metal
 
3
atom_style charge
 
4
dielectric      1.
 
5
 
 
6
variable type string "_charge"
 
7
 
 
8
boundary        s s f
 
9
# read in CNT
 
10
read_data  cnt_9_0_100${type}.init
 
11
lattice       diamond  3.6 # NOTE ???
 
12
pair_style  airebo 3.0
 
13
pair_coeff  * * ./CH.airebo C
 
14
mass          *  12.01
 
15
 
 
16
# PARAMETERS-----------------------------
 
17
variable dt equal 0.0005
 
18
variable L equal zhi-zlo
 
19
variable zhi equal zhi
 
20
variable zTip equal ${zhi}-2.0 # 2 4
 
21
variable zFree equal zhi
 
22
variable R equal 12.1/2
 
23
variable xhiFE equal 5.0*$R
 
24
variable xloFE equal -${xhiFE}
 
25
variable yhiFE equal $R
 
26
variable yloFE equal -${yhiFE}
 
27
variable zloFE equal zlo+10 # create fixed ghosts
 
28
variable zhiFE equal zhi+(zhi-${zloFE})/12*2
 
29
variable Lfree equal zhi-${zloFE}
 
30
variable nx equal 10 # 5
 
31
variable nz equal 14 # 12
 
32
print "Length $L [${zloFE}, ${zhiFE}] ${zTip}"
 
33
 
 
34
#variable E  equal 0.1  # bias 1.0
 
35
variable Vb equal 0.1 # 0.1 #0.5 #0.0 # bias
 
36
variable Vg equal 0.5 # 1.0 # 5.0 0.5 #50.0 # 0.5 # gate/modulation
 
37
print "bias voltage ${Vb}, gate voltage ${Vg}"
 
38
 
 
39
variable ng equal 20 # 80 # 10
 
40
variable nb equal 2 # 3
 
41
 
 
42
variable n equal 100000
 
43
variable s equal 250
 
44
# END -----------------------------------
 
45
 
 
46
region  TIP block INF INF INF INF ${zTip} INF  units box
 
47
group   TIP region TIP
 
48
 
 
49
#region  FIXED block INF INF INF INF INF ${zLoFE}  units box
 
50
#group   FIXED region FIXED
 
51
 
 
52
#group   FREE subtract all FIXED
 
53
 
 
54
region  feRegion block ${xloFE} ${xhiFE} ${yloFE} ${yhiFE} ${zloFE} ${zhiFE} units box
 
55
group internal region feRegion
 
56
group FIXED subtract all internal
 
57
fix FIX FIXED setforce 0 0 0
 
58
 
 
59
variable nAll   equal count(all)
 
60
variable nGhost equal count(all)-count(internal)
 
61
variable nTip   equal count(TIP)
 
62
print ">>> number of stationary ghosts:  ${nGhost} of ${nAll}"
 
63
print ">>> number of tip atoms :  ${nTip}"
 
64
 
 
65
neighbor        5. bin
 
66
neigh_modify    every 10 delay 0 check no
 
67
timestep        ${dt} 
 
68
thermo 100
 
69
variable tag string "electrostatic_bending"
 
70
 
 
71
# set charge on tip
 
72
variable C equal -0.025 
 
73
print "charge $C [e]"
 
74
variable c equal $C/${nTip}
 
75
set group TIP charge $c
 
76
 
 
77
# coupling
 
78
fix         AtC internal   atc electrostatic CNT_id.mat
 
79
fix_modify  AtC include atomic_charge
 
80
fix_modify  AtC internal_quadrature off
 
81
# note weights don't affect phi or f
 
82
fix_modify  AtC atom_weight constant internal 1.0
 
83
fix_modify  AtC extrinsic short_range off
 
84
#fix_modify  AtC atom_element_map eulerian 1
 
85
fix_modify  AtC control  momentum flux
 
86
 
 
87
fix_modify  AtC mesh create ${nx} 1 ${nz} feRegion f p f
 
88
 
 
89
# initial conditions
 
90
fix_modify  AtC  initial displacement x all 0.0
 
91
fix_modify  AtC  initial displacement y all 0.0
 
92
fix_modify  AtC  initial displacement z all 0.0
 
93
fix_modify  AtC  initial velocity x all 0.0
 
94
fix_modify  AtC  initial velocity y all 0.0
 
95
fix_modify  AtC  initial velocity z all 0.0
 
96
fix_modify  AtC  initial electric_potential all 0.0
 
97
 
 
98
# node sets
 
99
variable t equal 1.1*$R
 
100
fix_modify  AtC mesh create_nodeset tube -$t $t -$t $t ${zloFE} ${zFree} units box
 
101
fix_modify  AtC mesh create_nodeset lefttube -$t $t -$t $t ${zloFE} ${zloFE} units box
 
102
fix_modify  AtC mesh create_nodeset rbc INF INF INF INF ${zhiFE} ${zhiFE} units box
 
103
fix_modify  AtC mesh create_nodeset lbc INF INF INF INF ${zloFE} ${zloFE} units box
 
104
fix_modify  AtC mesh create_nodeset top ${xhiFE} ${xhiFE} INF INF INF INF units box
 
105
fix_modify  AtC mesh create_nodeset bot ${xloFE} ${xloFE} INF INF INF INF units box
 
106
 
 
107
# boundary conditions
 
108
fix_modify AtC  fix displacement x lbc 0.
 
109
fix_modify AtC  fix displacement y lbc 0.
 
110
fix_modify AtC  fix displacement z lbc 0.
 
111
fix_modify AtC  fix velocity     x lbc 0.
 
112
fix_modify AtC  fix velocity     y lbc 0.
 
113
fix_modify AtC  fix velocity     z lbc 0.
 
114
# ground
 
115
fix_modify AtC  fix electric_potential lbc 0
 
116
# bias
 
117
fix_modify AtC  fix electric_potential rbc ${Vb}
 
118
# gate
 
119
fix_modify AtC  fix electric_potential bot ${Vg} 
 
120
 
 
121
# run
 
122
compute CM TIP com
 
123
compute q all property/atom q
 
124
compute Q  all reduce sum c_q
 
125
compute FSUM all reduce sum fx fy fz
 
126
compute RSUM internal reduce sum fx fy fz
 
127
 
 
128
thermo_style    custom step etotal ke c_CM[1] c_CM[2] c_CM[3]  &
 
129
  c_Q  f_AtC[4]  f_AtC[5]  f_AtC[6] f_AtC[7] f_FIX[1] f_FIX[2] f_FIX[3] f_AtC c_FSUM[1] c_RSUM[1]
 
130
thermo          $s
 
131
fix_modify      AtC  output  ${tag}FE 100000000 full_text # $s full_text #binary
 
132
fix_modify      AtC  output  index step      
 
133
# NOTE not recognized as vector by paraview
 
134
variable uX atom x-f_AtC[1]
 
135
variable uY atom y-f_AtC[2]
 
136
variable uZ atom z-f_AtC[3]
 
137
variable rho atom mass*f_AtC[4]
 
138
dump  CONFIG all custom $s ${tag}.dmp id type x y z v_uX v_uY v_uZ v_rho
 
139
 
 
140
reset_timestep 0
 
141
log ${tag}.log
 
142
 
 
143
# [eV/A * A^2] --> [N m]
 
144
variable eV2J equal 1.60217646e-19
 
145
variable A2m  equal 1.e-10
 
146
thermo 10
 
147
timestep   0.0
 
148
min_modify      line quadratic
 
149
variable Vg equal 0.1
 
150
variable Lx equal 1.0
 
151
variable ng equal 4
 
152
#compute RSUM FREE reduce sum fx fy fz
 
153
#dump  CONFIG all custom 10000 ${tag}.dmp id type x y z c_U[1] c_U[2] c_U[3] fx fy fz
 
154
variable a equal 0
 
155
  variable i loop ${ng}
 
156
  label loop_i
 
157
    variable b equal ($i-1)*${Vg}/(${ng}-1)/${Lx}
 
158
 
 
159
    fix_modify AtC  fix electric_potential all linear 0 0 0 $b 0 $a 0
 
160
    min_style cg
 
161
    min_modify line quadratic 
 
162
    #minimize 0 0 100000 100000
 
163
    minimize 0 0 1000 1000
 
164
    min_style sd
 
165
    min_modify line backtrack 
 
166
    #minimize 0 0 100000 100000
 
167
    minimize 0 0 1000 1000
 
168
    fix_modify AtC output now 
 
169
 
 
170
    # u = F L^3 / 3 EI --> EI = F L^3 / 3 u
 
171
    variable u equal c_CM[1]
 
172
    variable uz equal c_CM[3]
 
173
#   variable F equal f_AtC[5]
 
174
#   variable Fz equal f_AtC[7]
 
175
    variable F equal c_RSUM[1]
 
176
    variable Fz equal c_RSUM[3]
 
177
    variable R equal $F-$C*$b
 
178
    variable Rz equal ${Fz}-$C*$a
 
179
    variable EI equal $F*${Lfree}*${Lfree}*${Lfree}/3./$u
 
180
    variable EI equal ${EI}*${eV2J}*${A2m}
 
181
    #print "flexural rigidity ${EI}  [Nm^2]  NOTE z force"
 
182
 
 
183
    print ">> V $b $a F $F ${Fz} u $u ${uz} c $c phi 0 EI ${EI} R $R ${Rz}"
 
184
    next i
 
185
  jump SELF loop_i
 
186