NWChem - Revision history

imported>Adam at 20:50, March 14, 2013

2013-03-14T20:50:27Z

← Older revision		Revision as of 20:50, March 14, 2013
Line 230:		Line 230:

	In addition to reducing the complexity of calculations in g(), using a lookup table also reduces the number of memory operation from 20 (x, y, z, expnt and rnorm for each i,j,k,l index) in the reference code to 12 (2 * 6 struct entries).		In addition to reducing the complexity of calculations in g(), using a lookup table also reduces the number of memory operation from 20 (x, y, z, expnt and rnorm for each i,j,k,l index) in the reference code to 12 (2 * 6 struct entries).

			===Symmetric matrices===
			Both the input, g_dens, and output, g_fock, matrices used in twoel() are symmetric matrices. This means that half of the memory accesses or updates to each are redundant and can be eliminated from the inner loops of twoel(). The required changes to the updates vary somewhat between the Coulomb and exchange contribution calculations. The Coulomb update contribution is performed by the following function:

			void update_coulomb(double g_value, int a, int b, int c, int d) {
			g_fock[a][b] += g_value * g_dens[c][d];
			}

			The symmetry in the Coulomb contributions is a subset of the symmetry in the g() function, and is therefore trivial to implement in the loops that have already been decomposed. Consider the updates in the full 8-way g() symmetry loop which has the restrictions i != j, k != l and i,j != k,l:

			update_coulomb(g_value, i, j, k, l);
			update_coulomb(g_value, i, j, l, k);
			update_coulomb(g_value, j, i, k, l);
			update_coulomb(g_value, j, i, l, k);
			update_coulomb(g_value, k, l, i, j);
			update_coulomb(g_value, k, l, j, i);
			update_coulomb(g_value, l, k, i, j);
			update_coulomb(g_value, l, k, j, i);

			The "a" and "b" variables are used as indexes into the g_fock output matrix, when they are swapped, the same value is accumulated into the corresponding symmetric entry. All updates which swap "a" and "b" can be removed as the loop has already been decomposed to insure that j > i and k > l, therefore the swapped values only update the lower triangle which is unnecessary. This leaves us with the following updates:

			update_coulomb(g_value, i, j, k, l);
			update_coulomb(g_value, i, j, l, k);
			update_coulomb(g_value, k, l, i, j);
			update_coulomb(g_value, k, l, j, i);

			Next we can use the symmetry of g_dens to remove the updates which swap the "c" and "d" arguments. Because this is using symmetry from the input matrix, we can't simply delete these updates. Instead, we can easily combine the symmetric updates into a single update. Recall that the Coulomb update is:
			g_fock[a][b] += g_value * g_dens[c][d];

			Instead of calling update_coulomb twice with swapped c and d arguments, the swap can be internalized to the function:
			g_fock[a][b] += g_value * g_dens[c][d] + g_value * g_dens[d][c];

			Because g_dens is symmetric, this can be substituted with:
			g_fock[a][b] += g_value * g_dens[c][d] + g_value * g_dens[c][d];

			Which is reduced to:
			g_fock[a][b] += 2.0 * g_value * g_dens[c][d];

			We can therefore create a new update function:

			void update_coulomb2(double g_value, int a, int b, int c, int d) {
			g_fock[a][b] += 2.0 * g_value * g_dens[c][d];
			}


			Using the new update_coulomb2() function leads to the final coulomb updates:

			update_coulomb2(g_value, i, j, k, l);
			update_coulomb2(g_value, k, l, i, j);

			The update_exchange() calls are more complicated to modify as their symmetry does not match the symmetry in g() as well. The exchange contributions are symmetric when "a" and "c" are swapped, or when "b" and "d" are swapped. This can only be easily used to get a two-fold reduction in updates which corresponds to the i,j != k,l loop restriction already in place. We did not fully take advantage of this symmetry to limit all updates to one triangle of the matrix. Instead, we eliminated all redundant updates but still allowed writing to both triangles. This required using a additional temporary matrix which was then post-processed to regenerate the correct symmetric values. The updates to the temporary matrix are performed in the following function:

			void update_exchange(double g_value, int a, int b, int c, int d) {
			t_fock[a][c] -= 0.50 * g_value * g_dens[b][d];
			}

			When not taking advantage of any symmetry, the following updates must be performed:

			update_exchange(g_value, i, j, k, l);
			update_exchange(g_value, i, j, l, k);
			update_exchange(g_value, j, i, k, l);
			update_exchange(g_value, j, i, l, k);
			update_exchange(g_value, k, l, i, j);
			update_exchange(g_value, k, l, j, i);
			update_exchange(g_value, l, k, i, j);
			update_exchange(g_value, l, k, j, i);

			Using the i,j = k,l symmetry of g(), we can eliminate half of the updates:

			update_exchange(g_value, i, j, k, l);
			update_exchange(g_value, i, j, l, k);
			update_exchange(g_value, j, i, k, l);
			update_exchange(g_value, j, i, l, k);

			After completing the N^4 loop, the temporary matrix is accumulated into the fock matrix as follows:
			for (i = 0; i < nfbn; i++)
			for (j = 0; j < nfbn; j++)
			g_fock[i][j] += t_fock[i][j] + t_fock[j][i];

			Although we did not implement it, it should be possible to further decompose the loops to ensure only one triangle of the matrix is updated.

	===Optimizations for large problem sizes===		===Optimizations for large problem sizes===

imported>Rishi at 19:18, March 8, 2013

2013-03-08T19:18:51Z

← Older revision		Revision as of 19:18, March 8, 2013
Line 261:		Line 261:
	====oneel()====		====oneel()====
	The calculation of h() in oneel() is unnecessary as the result never changes from one iteration to the next. Instead of recalculating h() to initialize the fock matrix each iteration, it can be calculated and stored once during initialization, then copied back into the fock matrix at the beginning of each iteration. The performance increase of this is small relative to the total execution time, but may still result noticeable differences for large problems.		The calculation of h() in oneel() is unnecessary as the result never changes from one iteration to the next. Instead of recalculating h() to initialize the fock matrix each iteration, it can be calculated and stored once during initialization, then copied back into the fock matrix at the beginning of each iteration. The performance increase of this is small relative to the total execution time, but may still result noticeable differences for large problems.

			== Parallelization ==
			This section is not complete yet. There is a good explanation in the [[media:BrandywineXStackReportQ2.pdf\|DynAX Q2 report]] here, but it is not yet incorporated into the wiki.
			=== Shared Memory ===
			==== OpenMP ====
			==== SWARM Single Node ====
			===Distributed Memory ===
			==== MPI + OpenMP ====
			==== SWARM Multiple Nodes ====
			.

imported>Rishi: /* SCF Module */

2013-03-07T02:58:33Z

SCF Module

← Older revision		Revision as of 02:58, March 7, 2013
Line 19:		Line 19:

	== SCF Module ==		== SCF Module ==
	[[Media:SCF.pptx \| Click here to view a Powerpoint Presentation on the SCF Module]]		[[Media:SCF.pptx \| Click here to view a Powerpoint Presentation on the SCF Module]] <p>
			[[media:SCF_text_document.docx\|Click here to view a Word document describing the SCF module]]

	== Description ==		== Description ==

imported>Rishi: /* Download */

2013-03-07T02:57:06Z

Download

← Older revision		Revision as of 02:57, March 7, 2013
Line 24:		Line 24:

	== Download ==		== Download ==
	~~Code will be available~~ here ~~shortly~~. It is licensed under the BSD license.		You can download the code [[media:Scf.tar.gz\|here]]. It is licensed under the BSD license.

	== Optimizations ==		== Optimizations ==

imported>Jfeo at 19:54, February 18, 2013

2013-02-18T19:54:36Z

← Older revision		Revision as of 19:54, February 18, 2013
Line 14:		Line 14:
	== Applications ==		== Applications ==

	There are many modules of [http://nwchem-sw.org NWChem]. In the DynAX XStack program, PNNL will extract ~~several~~ of these modules and produce simple, readable, serial modules. The DynAX team will produce parallel and optimized versions of these modeules and report on those optimizations.		There are many modules of [http://nwchem-sw.org NWChem]. In the DynAX XStack program, PNNL will extract two of these modules and produce simple, readable, serial modules. The DynAX team will produce parallel and optimized versions of these modeules and report on those optimizations.

	The first of these modules is the [http://www.nwchem-sw.org/index.php/Release61:Hartree-Fock_Theory_for_Molecules Self Consistent Field (SCF)].		The first of these modules is the [http://www.nwchem-sw.org/index.php/Release61:Hartree-Fock_Theory_for_Molecules Self Consistent Field (SCF)].

imported>Rishi at 19:00, February 7, 2013

2013-02-07T19:00:24Z

← Older revision		Revision as of 19:00, February 7, 2013
Line 14:		Line 14:
	== Applications ==		== Applications ==

	There are many modules of [http://nwchem-sw.org NWChem]. In the DynAX XStack program, PNNL will extract several of these modules and produce simple, readable, serial ~~kernels~~. The DynAX team will produce parallel and optimized versions of these ~~kernels~~ and report on those optimizations.		There are many modules of [http://nwchem-sw.org NWChem]. In the DynAX XStack program, PNNL will extract several of these modules and produce simple, readable, serial modules. The DynAX team will produce parallel and optimized versions of these modeules and report on those optimizations.

	The first of these ~~kernels~~ is the [http://www.nwchem-sw.org/index.php/Release61:Hartree-Fock_Theory_for_Molecules Self Consistent Field (SCF)].		The first of these modules is the [http://www.nwchem-sw.org/index.php/Release61:Hartree-Fock_Theory_for_Molecules Self Consistent Field (SCF)].

	== SCF ~~Kernel~~ ==		== SCF Module ==
	[[Media:SCF.pptx \| Click here to view a Powerpoint Presentation on the SCF ~~kernel~~]]		[[Media:SCF.pptx \| Click here to view a Powerpoint Presentation on the SCF Module]]

	== Description ==		== Description ==
Line 27:		Line 27:

	== Optimizations ==		== Optimizations ==
	[[Media:SCF_optimizations.pdf \| Presentation on the SCF ~~kernel~~ optimizations (PDF format)]]		[[Media:SCF_optimizations.pdf \| Presentation on the SCF module optimizations (PDF format)]]
	===Additional Work Skipping===		===Additional Work Skipping===

imported>Adam at 20:06, February 5, 2013

2013-02-05T20:06:14Z

@@ Line 237: / Line 237: @@
 All of the replacement functions are direct equivalents which only require minor changes to use (mostly additional arguments).
 ====oneel()====
 The calculation of h() in oneel() is unnecessary as the result never changes from one iteration to the next. Instead of recalculating h() to initialize the fock matrix each iteration, it can be calculated and stored once during initialization, then copied back into the fock matrix at the beginning of each iteration. The performance increase of this is small relative to the total execution time, but may still result noticeable differences for large problems.

imported>Rishi: /* Using a precomputed lookup table in g() */

2013-02-05T18:37:16Z

Using a precomputed lookup table in g()

← Older revision		Revision as of 18:37, February 5, 2013
Line 196:		Line 196:
	double ex; // exp(-facij * rab2)		double ex; // exp(-facij * rab2)
	double expnt; // expnt[i] + expnt[j]		double expnt; // expnt[i] + expnt[j]
	double x; // x[i] * expnt[i] + x[j] * expnt[j]) / ~~expntnt~~		double x; // x[i] * expnt[i] + x[j] * expnt[j]) / expnt
	double y; // y[i] * expnt[i] + y[j] * expnt[j]) / ~~expntnt~~		double y; // y[i] * expnt[i] + y[j] * expnt[j]) / expnt
	double z; // z[i] * expnt[i] + z[j] * expnt[j]) / ~~expntnt~~		double z; // z[i] * expnt[i] + z[j] * expnt[j]) / expnt
	double rnorm; // rnorm[i] * rnorm[j]		double rnorm; // rnorm[i] * rnorm[j]
	} g_lookup_table;		} g_lookup_table;

imported>Adam at 18:32, February 5, 2013

2013-02-05T18:32:06Z

← Older revision		Revision as of 18:32, February 5, 2013
Line 27:		Line 27:

	== Optimizations ==		== Optimizations ==
			[[Media:SCF_optimizations.pdf \| Presentation on the SCF kernel optimizations (PDF format)]]
	===Additional Work Skipping===		===Additional Work Skipping===

imported>Adam at 18:09, February 5, 2013

2013-02-05T18:09:46Z

← Older revision		Revision as of 18:09, February 5, 2013
Line 229:		Line 229:

	In addition to reducing the complexity of calculations in g(), using a lookup table also reduces the number of memory operation from 20 (x, y, z, expnt and rnorm for each i,j,k,l index) in the reference code to 12 (2 * 6 struct entries).		In addition to reducing the complexity of calculations in g(), using a lookup table also reduces the number of memory operation from 20 (x, y, z, expnt and rnorm for each i,j,k,l index) in the reference code to 12 (2 * 6 struct entries).

			===Optimizations for large problem sizes===
			====Matrix routines====
			The matrix multiplication and eigenvalue solver used in diagon() operates fast for small problem sizes, but as the problem size increases, diagon() can consume a larger percentage of the overall time. In fact, if the input problem is sparse enough, diagon() can scale worse than twoel() and consume the majority of the computation time.

			This situation can be easily improved by replacing the Dgemm() in integ.c with a call to a BLAS dgemm and by replacing rs() and rsg() in diagonalize.c with dsyge() and dsygv() from LAPACK.

			All of the replacement functions are direct equivalents which only require minor changes to use (mostly additional arguments).

			====oneel()====
			The calculation of h() in oneel() is unnecessary as the result never changes from one iteration to the next. Instead of recalculating h() to initialize the fock matrix each iteration, it can be calculated and stored once during initialization, then copied back into the fock matrix at the beginning of each iteration. The performance increase of this is small relative to the total execution time, but may still result noticeable differences for large problems.