Hide keyboard shortcuts

Hot-keys on this page

r m x p   toggle line displays

j k   next/prev highlighted chunk

0   (zero) top of page

1   (one) first highlighted chunk

1from sys import argv 

2from numpy import where, array, finfo, float64, arange, zeros, savetxt 

3from math import sin, cos, pi 

4from contextlib import suppress 

5 

6def strings(): 

7 """ 

8 Mango Strings XYZ file generation. 

9 

10 commandline parameters 

11 [radius]%[tolerance] [number of particles] [postion Phi] [position Theta] [magnetisation Phi] [magnetisation Theta] 

12 

13 Phi and Theta are the azimuthal angle and the polar angle for position and magnetisation, 

14 they are optional and default to aligning along the z axis. 

15 -h help 

16 -v verbose 

17 

18 """ 

19 if '-h' in argv: 

20 print(strings.__doc__) 

21 exit(0) 

22 

23 v = v_print("-v" in argv) 

24 with suppress(ValueError): 

25 del argv[argv.index('-v')] 

26 

27 radius = argv[1] 

28 if '%' in radius: 

29 r = radius.split('%') 

30 radius = float(r[0]) + float(r[0]) * float(r[1]) / 100 

31 r = f'{r[0]}%{r[1]}' 

32 else: 

33 radius = float(radius) 

34 r = radius 

35 

36 no_mol = int(argv[2]) 

37 no_mol_arr = arange(no_mol) 

38 

39 try: 

40 p, t = list(array([argv[3], argv[4]], dtype=float) * (pi / 180)) 

41 mp, mt = list(array([argv[5], argv[6]], dtype=float) * (pi / 180)) 

42 except IndexError: 

43 v.print('Using default direction (z)') 

44 p = mp = t = mt = 0 

45 

46 unit, munit = (r_unit(p, t), r_unit(mp, mt)) 

47 

48 v.print('Position Unit vector', unit) 

49 v.print('Magnetic Unit vector', munit) 

50 

51 start = - (unit * radius * (no_mol - 1)) 

52 

53 xyz = zeros((no_mol, 9)) 

54 xyz[:, 0:3] = start + (no_mol_arr[:, None] * radius * 2 * unit) 

55 xyz[:, 6:9] += munit 

56 

57 v.print(xyz) 

58 

59 header = f' {no_mol}\n # Generated with mango strings generator' 

60 savetxt(f'string_{no_mol}r{r}.xyz', xyz, fmt='MNP {}'.format(9 * ' % .8e'), comments='', header=header) 

61 

62 

63def r_unit(p, t): 

64 unit = array([sin(t) * cos(p), sin(t) * sin(p), cos(t)]) 

65 return where(unit < finfo(float64).eps, where(unit > -finfo(float64).eps, 0, unit), unit) 

66 

67 

68class v_print(): 

69 

70 def __init__(self, verbose): 

71 if verbose: 

72 self.print = self._print 

73 else: 

74 self.print = self._dont_print 

75 

76 @staticmethod 

77 def _print(*a, **k): 

78 return print(*a, **k) 

79 

80 @staticmethod 

81 def _dont_print(*args): 

82 pass 

83 

84 

85if __name__ == '__main__': 

86 strings()