WRFSource Code Changes

What I Needed and Why It Is a Little Complicated

Earlier, as part of some work to assess wind for potential wind farm sites, I found that I wanted to use hourly averages for wind, not instantaneous such as WRF produces. And I need the wind speeds to be at specific heights above ground level (AGL) because wind turbines are. So I had WRF output very often and processed the huge data files to interpolate to the heights I wanted and calculate averages. Now that I am doing research for wind forecasting, I again want average wind speed over a time period; this time I want 10 minute averages to use to compare to the 10 minute average wind speeds reported by a sodar. Real winds are messy and the sodar's averaging helps even things out, but it seems appropriate to use 10 minute WRF averages if I am using 10 minute sodar averages to verify the results. So I needed a better solution than huge data files that require a lot of post processing. And I made one. I changed WRF so that it processes the data before it outputs it.

For creating averages, two steps are required: sum the values and divide by the count of the values. Not so fast. The values need to be comparable to each other. In this case, that means two things: they should be at the same height and they should be for equal length time periods. WRF uses pressure values for its vertical coordinate. Pressure changes over time and location, so converting to height AGL has to be done before summing. That needs to be done within WRF unless all the pressure data is also accumulated and reported.  So the code I added to WRF has to interpolate to the desired heights before summing the horizontal components, u (west-east) and v (south-north). So I need two new arrays in WRF, one for each component. Each array has to have three dimensions: west-east dimension, south-north dimension, and a vertical plane for each of the desired heights for which I want an average. I want to be flexible enough to have multiple heights. Only a single variable is needed for the count of the summed values because all three components in all three dimensions have the same number of values accumulated in the sums.

All the summed values should be for equal length time steps. Alternatively, the values can be weighted based on the relative length of the time step. I don't do this. Not yet. But no worries; I am currently using fixed length time steps. When I do more with adaptive time steps, I should make changes to weight things accordingly. The adapted time steps won't be vastly different from each other, but I might as well do things right; this is suppose to be science!

But I lied. I actually have four new arrays added to WRF. The third is for the vertical wind speed component. I can use this as data to feed into Fluent (see Initialize Fluent With WRF Data), but it is not used otherwise because I don't think wind turbines can capture energy from the vertical component of the wind and besides, it is usually a quite small value. The fourth array is for summing the combined horizontal wind speed before averaging. An average horizontal wind speed can not be calculated from the average of the U and V components. This example shows why: Suppose we have two data points. The first wind is from the NW at 7 m/s, so its U component is +5 m/s and the V is -5 m/s. The second wind is from the SW also at 7 m/s, so the U is also +5 m/s and the V is +5 m/s. We can use the average U, 5, and the average V, zero, to say that the average direction is from the west (u=5, v=0). That works. But if we combine the average U and the average V to get an average horizontal speed we get a 5 m/s speed from two 7 m/s winds. That does not work. We have to calculate an average horizontal speed before we accumulate a sum for it, separately from the components.

It would be nice to have WRF output the actual averages, but to do so would require extra processing and double the number of arrays, at least it would when I get to adaptive time steps. I want to keep my code changes to a minimum because WRF has new versions and bug fix releases in which I have to remake those changes. So it worked out best to let WRF use its usual routines to output the arrays of summed values and the count of those sums. It is a minor postprocessing step to do the final division of summed values by count of those values. I just have to reset the summing arrays and the count to zero after each time they are output. All my code changes go near the end of the one source file "solve_em.F" in the "dyn_em" directory.

Some Implementation Details

To create the four new arrays mentioned above, I make changes to the registry and add them as "state" variables. The horizontal dimensions are the same as the WRF grid, so that part was easy to specify. But the vertical is different. I didn't want to hard code the number of heights that can be averaged and there is no reason to tie it to the number of vertical levels WRF is using in its grid, so I created another dimension. This can be done with a "dimspec" variable in the registry. Actually, with newer versions of WRF, the dimspecs are defined in a file called "registry.dimspec," which gets included into the registry. I create a dimspec variable called "max_sums" that defines the number of heights. Then by specifying max_sums in the domains portion of the namelist.input file, the size of that dimension is set. The actual heights for which to calculate averages are also specified in the domains section of the namelist with a new variable called "fixed_sum_zs."

It might just be me, but I had trouble getting access to a couple variables I needed. One of those was max_sums. So I created another variable that I could use called "use_sums," which is specified in the domains section of the namelist. It tells how many of the fixed_sum_zs values to use when interpolating to heights and accumulating sums. For flexibility, I allowed a different number to be specified for each domain, but I'm not sure why I would want to do so other than to say zero for heights for outer domains. The other new variable defines how often to reset the summing array to zero, so of course it needs to match the output frequency of the history file. I couldn't find easy access to that information, so I had to create "sum_reset," which is also specified in the domains portion of the namelist. Since history output can differ per domain, I allow sum_reset to have different values for each domain. For my work, I want wind information output at 10 minute intervals, which is more frequent than I need other WRF information output, so I have the wind averaging information output into an auxiliary history file, number 6.

In case you didn't notice, let me point out that all this wind averaging stuff is run time configurable. There is no need to make source code changes and recompile if I want to average at a different set of heights.

Here are the changes that I made to the Registry files:

And here are the namelist.input entries I use to get averages for the same 10 heights at which my sodar reports (I only want averaging for my innermost domain, #3):

I don't really want to give out the actual source code I created to interpolate to the heights, sum the values, and reset the arrays, but if you want to write your own code, all the above information should give you a good running start.