ViewVC Help

View File | Revision Log | Show Annotations | View Changeset | Root Listing

root/OpenMD/trunk/src/applications/hydrodynamics/AnalyticalModel.cpp

Comparing trunk/src/applications/hydrodynamics/AnalyticalModel.cpp (file contents):
Revision 906 by tim, Fri Mar 17 23:20:35 2006 UTC vs.
Revision 1879 by gezelter, Sun Jun 16 15:15:42 2013 UTC

#	Line 6 | Line 6
6		* redistribute this software in source and binary code form, provided
7		* that the following conditions are met:
8		*
9	<	* 1. Acknowledgement of the program authors must be made in any
10	<	*    publication of scientific results based in part on use of the
11	<	*    program.  An acceptable form of acknowledgement is citation of
12	<	*    the article in which the program was described (Matthew
13	<	*    A. Meineke, Charles F. Vardeman II, Teng Lin, Christopher
14	<	*    J. Fennell and J. Daniel Gezelter, "OOPSE: An Object-Oriented
15	<	*    Parallel Simulation Engine for Molecular Dynamics,"
16	<	*    J. Comput. Chem. 26, pp. 252-271 (2005))
17	<	*
18	<	* 2. Redistributions of source code must retain the above copyright
9	>	* 1. Redistributions of source code must retain the above copyright
10		*    notice, this list of conditions and the following disclaimer.
11		*
12	<	* 3. Redistributions in binary form must reproduce the above copyright
12	>	* 2. Redistributions in binary form must reproduce the above copyright
13		*    notice, this list of conditions and the following disclaimer in the
14		*    documentation and/or other materials provided with the
15		*    distribution.
#	Line 37 | Line 28
28		* arising out of the use of or inability to use software, even if the
29		* University of Notre Dame has been advised of the possibility of
30		* such damages.
31	+	*
32	+	* SUPPORT OPEN SCIENCE!  If you use OpenMD or its source code in your
33	+	* research, please cite the appropriate papers when you publish your
34	+	* work.  Good starting points are:
35	+	*
36	+	* [1]  Meineke, et al., J. Comp. Chem. 26, 252-271 (2005).
37	+	* [2]  Fennell & Gezelter, J. Chem. Phys. 124, 234104 (2006).
38	+	* [3]  Sun, Lin & Gezelter, J. Chem. Phys. 128, 234107 (2008).
39	+	* [4]  Kuang & Gezelter,  J. Chem. Phys. 133, 164101 (2010).
40	+	* [5]  Vardeman, Stocker & Gezelter, J. Chem. Theory Comput. 7, 834 (2011).
41		*/
42		#include "applications/hydrodynamics/AnalyticalModel.hpp"
43	<	#include "applications/hydrodynamics/Spheric.hpp"
44	<	#include "applications/hydrodynamics/Ellipsoid.hpp"
43	>	#include "hydrodynamics/Sphere.hpp"
44	>	#include "hydrodynamics/Ellipsoid.hpp"
45		#include "applications/hydrodynamics/CompositeShape.hpp"
46		#include "math/LU.hpp"
47	<	namespace oopse {
48	<	bool AnalyticalModel::calcHydroProps(Spheric* spheric, double viscosity, double temperature) {
49	<
49	<	double radius = spheric->getRadius();
50	<	HydroProps props;
51	<	props.center =V3Zero;
52	<	double Xitt  = 6.0 * NumericConstant::PI * viscosity * radius;
53	<	double Xirr = 8.0 * NumericConstant::PI * viscosity * radius * radius * radius;
54	<	props.Xi(0, 0) = Xitt;
55	<	props.Xi(1, 1) = Xitt;
56	<	props.Xi(2, 2) = Xitt;
57	<	props.Xi(3, 3) = Xirr;
58	<	props.Xi(4, 4) = Xirr;
59	<	props.Xi(5, 5) = Xirr;
47	>	namespace OpenMD {
48	>
49	>	bool AnalyticalModel::calcHydroProps(Shape* shape, RealType viscosity, RealType temperature) {
50
51	<	const double convertConstant = 6.023; //convert poise.angstrom to amu/fs
52	<	props.Xi *= convertConstant;
53	<	Mat6x6d XiCopy = props.Xi;
54	<	invertMatrix(XiCopy, props.D);
55	<	double kt = OOPSEConstant::kB * temperature;
56	<	props.D *= kt;
57	<	props.Xi *= OOPSEConstant::kb * temperature;
58	<
59	<	setCR(props);
60	<	setCD(props);
61	<
62	<	return true;
51	>	HydroProp* props;
52	>	Sphere* sphere = dynamic_cast<Sphere*>(shape);
53	>	if (sphere != NULL) {
54	>	props = sphere->getHydroProp(viscosity, temperature);
55	>	setCR(props);
56	>	setCD(props);
57	>	return true;
58	>	} else {
59	>	Ellipsoid* ellipsoid = dynamic_cast<Ellipsoid*>(shape);
60	>	if (ellipsoid != NULL) {
61	>	props = ellipsoid->getHydroProp(viscosity, temperature);
62	>	setCR(props);
63	>	setCD(props);
64	>	return true;
65	>	} else {
66	>	CompositeShape* composite = dynamic_cast<CompositeShape*>(shape);
67	>	if (composite != NULL) {
68	>	//           props = composite->getHydroProp(viscosity, temperature);
69	>	//           setCR(props);
70	>	//           setCD(props);
71	>	//           return true;
72	>	return false;
73	>	} else {
74	>	sprintf( painCave.errMsg,
75	>	"Could not figure out what kind of shape this is!\n");
76	>	painCave.severity = OPENMD_ERROR;
77	>	painCave.isFatal = 1;
78	>	simError();
79	>	return false;
80	>	}
81	>	}
82	>	}
83	>	}
84
85	+	void AnalyticalModel::writeBeads(std::ostream& os) {
86	+	os << "1\n";
87	+	os << "Generated by Hydro\n";
88	+	Vector3d pos = sd_->getPos();
89	+	os << sd_->getType() << "\t" << pos[0] << "\t" << pos[1] << "\t" << pos[2] << std::endl;
90	+	}
91		}
75	–
76	–	/**
77	–	* calculate the ratio of friction coeffiction constant between ellipsoid and spheric
78	–	* with same volume.
79	–	* @param m
80	–	* @param n
81	–	* @note
82	–	* Reference:
83	–	*
84	–	* (1) Victor A. Bloomfield, On-Line Biophysics Textbook, Volume: Separations and Hydrodynamics
85	–	* Chapter 1,Survey of Biomolecular Hydrodynamics
86	–	* http://www.biophysics.org/education/vbloomfield.pdf
87	–	* (2) F. Perrin , J. Phys. Radium, [7] 5, 497-511, 1934
88	–	* (3) F. Perrin, J. Phys. Radium, [7] 7, 1-11, 1936
89	–	*/
90	–	bool AnalyticalModel::calcHydroProps(Ellipsoid* ellipsoid, double viscosity, double temperature) {
91	–	double ft;
92	–	double fra;
93	–	double frb;
94	–	double a = ellipsoid->getA();
95	–	double b = ellipsoid->getB();
96	–	double q = a/b; //?
97	–	if (q > 1.0) {//prolate
98	–	ft = sqrt(1-q*q)/(pow(q, 2.0/3.0)*log((1 + sqrt(1-q*q))/q));
99	–	fra = 4*(1-q*q)/(3*(2 - 2*pow(q, 4.0/3.0)/ft)); //not sure
100	–	frb = 4*(1-q*q*q*q) /(3*q*q*(2*pow(q, -2.0/3.0)*(2-q*q)/ft-2));
101	–	} else {//oblate
102	–	ft = sqrt(1-q*q)/(pow(q, 2.0/3.0)*atan(sqrt(q*q-1)));
103	–	fra = 4*(1-q*q)/(3*(2 - 2*pow(q, 4.0/3.0)/ft)); //not sure
104	–	frb = 4*(1-q*q*q*q) /(3*q*q*(2*pow(q, -2.0/3.0)*(2-q*q)/ft-2));
105	–	}
106	–
107	–	double radius = pow(a*a*b, 1.0/3.0);
108	–	HydroProps props;
109	–	double Xitt  = 6.0 * NumericConstant::PI * viscosity * radius;
110	–	double Xirr = 8.0 * NumericConstant::PI * viscosity * radius * radius * radius;
111	–	props.Xi(0, 0) = Xitt;
112	–	props.Xi(1, 1) = Xitt;
113	–	props.Xi(2, 2) = Xitt;
114	–	props.Xi(3, 3) = Xirr;
115	–	props.Xi(4, 4) = Xirr;
116	–	props.Xi(5, 5) = Xirr;
117	–
118	–	const double convertConstant = 6.023; //convert poise.angstrom to amu/fs
119	–	props.Xi *= convertConstant;
120	–	props.Xi(0,0) *= ft;
121	–	props.Xi(1,1) *= ft;
122	–	props.Xi(2,2) *= ft;
123	–	props.Xi(3,3) *= fra;
124	–	props.Xi(4,4) *= fra;
125	–	props.Xi(5,5) *= frb;
126	–
127	–	Mat6x6d XiCopy = props.Xi;
128	–	XiCopy /= OOPSEConstant::kb * temperature;
129	–	invertMatrix(XiCopy, props.D);
130	–	double kt = OOPSEConstant::kB * temperature;
131	–	props.D *= kt;
132	–
133	–	setCR(props);
134	–	setCD(props);
135	–
136	–	return true;
137	–	}
138	–
139	–	bool AnalyticalModel::calcHydroProps(CompositeShape* compositexShape, double viscosity, double temperature) {
140	–	return false;
141	–	}
142	–
143	–
144	–
145	–	}

Comparing trunk/src/applications/hydrodynamics/AnalyticalModel.cpp (property svn:keywords):
Revision 906 by tim, Fri Mar 17 23:20:35 2006 UTC vs.
Revision 1879 by gezelter, Sun Jun 16 15:15:42 2013 UTC

#	Line 0 | Line 1
1	+	Author Id Revision Date

Diff Legend

–	Removed lines
+	Added lines
<	Changed lines
>	Changed lines