text:
# Script to read IANA time zone files and create C header files defining the rules and time zones
module FileProcessor # Clears definitions
# The following is necessary to overload the operators for the new datatypes
using Base
using Nullables
using Dates
# Establish the directory containing the time zone files and the directory to store the
#  header file
# Open preferentially on the Z: drive, which will be the case when working at home
if isdir("Z:")
  TZFileDir = "Z:\\Ultimate Datetime Datatype\\TZ Database\\TZ Database\\"
  TZHeaderDir = "Z:\\Ultimate Datetime Datatype\\Ultimate Datetime\\LocalDatetimeDLL\\"
# Check for local directory on Windows machine
elseif isdir("C:")
  TZFileDir = "C:\\Users\\jay\\Documents\\Ultimate Datetime Datatype\\TZ Database\\TZ Database\\"
  TZHeaderDir = "C:\\Users\\jay\\Documents\\Ultimate Datetime Datatype\\Ultimate Datetime\\LocalDatetimeDLL\\"
# Check for OS / X
elseif isdir("/Volumes")
  TZFileDir = "/Volumes/Shared/Ultimate Datetime Datatype/TZ Database/TZ Database/"
  TZHeaderDir = "/Volumes/Shared/Ultimate Datetime Datatype/Ultimate Datetime/LocalDatetimeDLL/"
end
# Create a helper function to add a string to the end of a string array in a dictionary
function addString(nameDict::Dict{String,Array{String}}, key::AbstractString, incValue::AbstractString)
  if (haskey(nameDict, key))
    push!(nameDict[key], incValue)
  else
    nameDict[key] = [incValue]
  end
end
#
# Establish an array of TZ file names
TZFile = ["africa", "antarctica", "asia", "australasia", "europe", "northamerica", "southamerica"]
# Create a rule structure
struct TZRule
  name::String
  fromYear::UInt32 # First year rule is applicable
  toYear::UInt32 # Last year rule is applicable
  transMonth::UInt8 # Month the transition to or from daylight savings time takes place
  afterDayOfMonth::UInt8 # The day of the month the transition takes place, if week == 0,
                         #  or the day on which to begin counting if weekNumber > 0 &&
                         #  weekNumber < 6.
  weekNumber::UInt8      # 1 = first week of the month, 2 = second ..., 6 = last week of the month
  dayOfWeek::UInt8       # 0 = Sunday, 1 = Monday, etc.
  hourOfTransition::UInt8 # Hour of the day the transition is to take place
  minuteOfTransition::UInt8 # Minute of the day the transition is to take place
  secondOfTransition::UInt8 # Second of the day the transition is to take place
  hundredthOfTransition::UInt8 # Hundredth of a second of the day the transition is to take place
  basisOfTransition::UInt8 # 0 = universal time, 1 = local standard time, 2 = local wall time
  hoursSaved::UInt8 # Hours saved vs standard time after the transition
  minutesSaved::UInt8 # Minutes saved vs standard time after the transition
  abbrevChars::String # Characters used to modify the time zone abbreviation when the rule is in effect
end
# Create a time zone structure
struct TimeZone
  name::String
  # GMT Offset is the amount of time to add to UTC to get the standard time in this time zone
  isNegativeGMTOffset::UInt8 # 0 = positive offset, 1 = negative offset
  hoursOfGMTOffset::UInt8 # Number of hours of GMT offset
  minutesOfGMTOffset::UInt8 # Number of minutes of GMT offset
  secondsOfGMTOffset::UInt8 # Number of seconds of GMT offset
  hundredthsOfGMTOffset::UInt8 # Number of hundredths of a second of GMT Offset
  # A time zone can either have no daylight savings time rule, 1 daylight savings time rule, or
  #  a specified time offset.
  # If the rule filed is blank, use the offset.  When no rule is specified, set the offset to zero,
  #  to handle the second case, above.
  daylightSavingsTimeRule::String
  hoursSaved::UInt8 # Number of hours saved vs standard time after the transition
  minutesSaved::UInt8 # Number of minutes saved vs standard time after the transition
  # The following field contains the time zone abbreviation.  If it includes a slash,
  #  the field contains the abbreviations for both standard and daylight savings time.
  # If it includes %s, this is a place holder for the characters provided by the daylight
  #  savings time rules.
  abbrevChars::String
  # The until datetime indicates the date and time until which the time zone is effective,
  #  using the GMT offset and daylight savings time rules of the time zone.
  #  The default basis of this datetime is local wall time.
  # The basis can be modified to local standard time or to universal time.
  untilYear::UInt32 # Year of the datetime until which the time zone is effective.  No need to consider gigayear.
  untilMonth::UInt8 # Month of the datetime until which the time zone is effective.
  untilDayOfMonth::UInt8 # Day of the datetime until which the time zone is effective.
  untilHour::UInt8 # Hour of the datetime until which the time zone is effective.
  untilMinute::UInt8 # Minute of the datetime until which the time zone is effective.
  untilSecond::UInt8 # Second of the datetime until which the time zone is effective.
  untilHundredth::UInt8 # Hundredth of a second of the datetime until which the time zone is effecttive.
  untilBasis::UInt8 # 0 = universal time, 1 = local standard time, 2 = local wall time
end
# Read each of the TZ files and create the rule header file
# Declare the array to contain the rules
rulesArray = Array{String}(undef,3000,9)
nrules = 0
maxRNameLength = 0
maxAbbrevLength = 0
# The following regular expression will be used to identify candidate
#  fractional second GMT Offset specifications withing the comments in
#  the files.  The IANA software does not support fractional seconds, so
#  the values were rounded to whole seconds.
#
# We could manually edit the file to make the changes, but this would
#  have to be done every time we received updated files.
#
# Instead, we create an array of candidate fractional specifications
#  and substitute them if a matching (rounded) GMTOff is located
#
fracSecSpec = r"-?[0-9][0-9]?:[0-9][0-9]?:[0-9][0-9]?\.[0-9][0-9]?"
candidateFracSpec = Array{Int}(undef, 1000, 8)
numCandidates = 0
numMatches = 0
for file in TZFile
  infile = open(TZFileDir * file,"r")
  lines = readlines(infile)
  for line in lines
    # Search for fractional second specification on comment lines
    fracsec = match(fracSecSpec, line)
    if (length(line) > 1 && line[1] == '#' && fracsec != nothing)
      # Process as if it is a GMTOff field
      global numCandidates += 1
      # Check for a negative sign
      if (fracsec.match[1] == '-')
        isneg = 1
        gmtoff = fracsec.match[2:end]
      else
        isneg = 0
        gmtoff = fracsec.match
      end
      # Default all values to 0
      hour = 0
      minute = 0
      second = 0
      hundredth = 0
      # A '-' indicates 0
      if (gmtoff[1] != '-')
        x = split(gmtoff,":")
        hour = parse(Int,x[1])
        if (length(x) > 1)
          minute = parse(Int,x[2])
          if (length(x) == 3)
            # Split again to locate fractional seconds
            y = split(x[3],".")
            second = parse(Int,y[1])
            if (length(y) == 2)
              hundredth = parse(Int,y[2])
              if (hundredth < 10)
                hundredth *= 10
              end
            end
          end
        end
      end
      # Round the value to the next higher second.  Do this instead of
      #  actual rounding, since the rounding is not done consistently in the file
      rndhour = hour
      rndminute = minute
      rndsecond = second
      if (hundredth >= 0)  # Had been: if (hundredth >= 50)
        rndsecond += 1
        if (rndsecond == 60)
          rndsecond = 0
          rndminute += 1
          if (rndminute == 60)
            rndminute = 0
            rndhour += 1
          end
        end
      end
      # Store the values in the candidate fractional second GMTOFF array
      candidateFracSpec[numCandidates,1] = 0 # Number of times this entry has been matched
      candidateFracSpec[numCandidates,2] = rndhour
      candidateFracSpec[numCandidates,3] = rndminute
      candidateFracSpec[numCandidates,4] = rndsecond
      candidateFracSpec[numCandidates,5] = hour
      candidateFracSpec[numCandidates,6] = minute
      candidateFracSpec[numCandidates,7] = second
      candidateFracSpec[numCandidates,8] = hundredth
    end
    # The following line defends against Unicode characters in the string, which can generate StringIndexError
    if (length(line) > 4 && line[1:1] == "R" && line[2:4] == "ule")
      global nrules += 1
      nsubs = 0
      for sub in split(line)
        nsubs += 1
        # Check for maximum name length
        if (nsubs == 2)
          global maxRNameLength = max(maxRNameLength, length(sub))
        end
        # Check for maximum abberviation characters length
        if (nsubs == 10)
          global maxAbbrevLength = max(maxAbbrevLength, length(sub))
        end
        # Sometimes a comment and another word or phrase is added at the end of a rule
        if (sub == "#")
          break
        elseif (nsubs > 1)
          rulesArray[nrules, nsubs - 1] = sub
        end
      end
    end
  end
  close(infile)
end
print("Processed " * string(nrules) * " rule versions.\n")
print("Longest rule name: " * string(maxRNameLength) * " characters.\n")
print("Longest abberviation modifier: " * string(maxAbbrevLength) * " characters.\n")
print("Number of candidates for fractional second GMTOff specification: " *
          string(numCandidates) * "\n")
#
# Create an array of tzRuleVersions from the rulesArray
# Create a tzRules dictionary, which indexes into the tzRuleVersions array
#
tzRuleVersions = Array{TZRule}(undef,nrules)
tzRules = Dict{String,Array{String}}()
# Create an array of month abbreviations
months = ["Jan", "Feb", "Mar", "Apr", "May", "Jun", "Jul", "Aug", "Sep", "Oct", "Nov", "Dec"]
# Createt an array of day abbreviations
days = ["Sun", "Mon", "Tue", "Wed", "Thu", "Fri", "Sat"]
daysj = ["Mon", "Tue", "Wed", "Thu", "Fri", "Sat", "Sun"]
# Process the rules array
maxSavedHours = 0
maxSavedMinutes = 0
maxTZ = ""
previousName = ""
for i = 1:nrules
  addString(tzRules,rulesArray[i,1], string(i-1))
  # **** FROM field ****
  # Check whether any of the from values are not integers between 1900 and 2100
  from = 0
  v = tryparse(Int,rulesArray[i,2])
  if (v == nothing)
    print("From field not an integer - " * rulesArray[i,2])
  else
    v = parse(Int,rulesArray[i,2])
    if (v < 1900 || v > 2100)
      print("From field out of range: " * string(v))
    else
      from = v
    end
  end
  # **** TO field ****
  # Check whether any of the to values are not integers between 1900 and 2100
  to = 0
  v = tryparse(Int,rulesArray[i,3])
  if (v == nothing)
    #Check for 'only'
    if (rulesArray[i,3] == "only")
      to = from
    elseif (rulesArray[i,3] == "max")
      to = 1000000000 # Will be translated to 99,999,999,999
    else
      print("To field is not valid - " * rulesArray[i,3] * "\n")
    end
  else
    v = parse(Int,rulesArray[i,3])
    if (v < 1900 || v > 2100)
      print("To field out of range: " * string(v) * "\n")
    else
      to = v
    end
  end
  # **** Ignore the TYPE field, as it is blank in all cases
  if (rulesArray[i,4] != "-")
    print(rulesArray[i] * "\n")
  end
  # **** IN field
  month = 0
  for j in eachindex(months)
    if (rulesArray[i,5] == months[j])
      month = j
      break
    end
  end
  # Check for invalid month
  if (month == 0)
    print("Invalid month - " * rulesArray[i,5] * "\n")
  end
  # **** ON field
  dom = 0
  week = 0
  day = 7
  v = tryparse(Int,rulesArray[i,6])
  if (v == nothing)
    if (occursin(">=",rulesArray[i,6]))
      # day>=dom specification
      for j in eachindex(days)
        if (rulesArray[i,6][1:3] == days[j])
          day = j - 1
          break
        end
      end
      # Check for invalid day
      if (day == 7)
        print("Invalid day - " * rulesArray[i,6][1:3] * "\n")
      end
      # Day of month begins in character 6
      dom = parse(Int, rulesArray[i,6][6:end])
      # Handle >= 1, 8, 15 and 22 as 1st, 2nd, 3rd, 4th, respectively to account for skipped days,
      #  even though it is unlikely there was a rule processed during a transition month
      if (dom % 7 == 1)
        week = div(dom,7) + 1
        dom = dom % 7
      end
    elseif (occursin("last",rulesArray[i,6]))
      # lastday specification
      week = 6
      for j in eachindex(days)
        if (rulesArray[i,6][5:7] == days[j])
          day = j - 1
          break
        end
      end
      # Check for invalid day
      if (day == 7)
        print("Invalid day - " * rulesArray[i,6][1:3] * "\n")
      end
    else
      print("Invalid day specification: " * rulesArray[i,6])
    end
  else
    # Day of month specification
    dom = parse(Int,rulesArray[i,6])
  end
  # **** AT field
  # Default to midnight
  hour = 0
  minute = 0
  second = 0
  hundredth = 0
  # Determin the AT time basis by checking the last character of the AT field
  basis = 2
  at = rulesArray[i,7]
  v = tryparse(Int,string(rulesArray[i,7][end]))
  if (v == nothing)
    if (rulesArray[i,7[end]] == 'u')
      basis = 0
    elseif (rulesArray[i,7][end] == 's')
      basis = 1
    else
      basis = 2
    end
    # Remove the last character for subsequent processing
    at = rulesArray[i,7][1:end-1]
  else
    at = rulesArray[i,7]
  end
  # A - indicates midnight
  if (at[1] != '-')
    x = split(at,":")
    hour = parse(Int,x[1])
    if (length(x) > 1)
      minute = parse(Int,x[2])
      if (length(x) == 3)
        # Split again to locate fractional seconds
        y = split(x[3],".")
        second = parse(Int,y[1])
        if (length(y) == 2)
          hundredth = parse(Int,y[2])
          if (hundredth < 10)
            hundredth *= 10
          end
        end
      end
    end
  end
  # **** SAVE field
  hourssaved = 0
  minutessaved = 0
  save = rulesArray[i,8]
  x = split(save,":")
  hourssaved = parse(Int,x[1])
  if (length(x) == 2)
    minutessaved = parse(Int,x[2])
  end
  if (60*hourssaved + minutessaved > 60*maxSavedHours + maxSavedMinutes)
    global maxSavedHours = hourssaved
    global maxSavedMinutes = minutessaved
    global maxTZ = rulesArray[i,1]
  end
  # **** LETTER field
  letters = (rulesArray[i,9] != "-" ? rulesArray[i,9] : "")
  # Create a TZRule
  tzRuleVersions[i] = TZRule(rulesArray[i,1], from, to, month, dom, week, day, hour, minute,
                  second, hundredth, basis, hourssaved, minutessaved, letters)
end
# Create another rule versions array, which has all of the versions of a rule contiguous
tzRuleVersions2 = Array{TZRule}(undef,nrules)
v2ind = 0
for key in keys(tzRules)
    val = tzRules[key]
    for j = 1:length(val)
        vind = parse(Int, val[j]) + 1
        global v2ind += 1
        tzRuleVersions2[v2ind] = tzRuleVersions[vind]
    end
end
# Reconstruct tzRules
tzRules = Dict{String,Array{String}}()
for i = 1:nrules
    addString(tzRules, tzRuleVersions2[i].name, string(i-1))
end
# Create a dictionary containing the rule sets
#   Contains the indices of the versions of a rule that are applicable for a time range
struct AppRuleVersionSet
    startYear:: UInt32
    endYear:: UInt32
    # Assumes up to 4 rule versions can be applicable at any time
    versionIndexArray:: Array{UInt32}
end
# Create a helper function to add a set of rule versions to the end of a rule version set array in a dictionary
function addRuleSet(ruleSetDict::Dict{String,Array{AppRuleVersionSet}}, key::AbstractString, incValue::AppRuleVersionSet)
  if (haskey(ruleSetDict, key))
    push!(ruleSetDict[key], incValue)
  else
    ruleSetDict[key] = [incValue]
  end
end
#
ruleVersionSet = Dict{String,Array{AppRuleVersionSet}}()
# Loop through the rule versions deriving the applicable rule sets
#  tzRules contains the indices of the rules
for key in keys(tzRules)
    firstInd =  parse(Int,tzRules[key][1]) + 1
    lastInd = parse(Int,tzRules[key][end]) + 1
    firstYear = 1000000000
    lastYear = 0
    # Find the earliest starting year and the latest end year for any rule version
    for i = firstInd:lastInd
        if (tzRuleVersions2[i].fromYear < firstYear)
            firstYear = tzRuleVersions2[i].fromYear
        end
        if (tzRuleVersions2[i].toYear > lastYear && tzRuleVersions2[i].toYear != 1000000000)
            lastYear = tzRuleVersions2[i].toYear
        end
    end
    # Now find all of the rule versions that start in the same year
    indexArray = Array{Int}(undef,0)
    for i = firstInd:lastInd
        if (tzRuleVersions2[i].fromYear == firstYear)
            # Change from 1 based array to 0 based Array
            push!(indexArray,i-1)
        end
    end
    # Loop through all of the years from the first to the last year looking for
    #   changes in applicable rules
    begYear = firstYear
    for i = firstYear+1:lastYear+1
        upperLocated = false
        deletedIndices = Array{Int}(undef,0)
        # Check for rule versions beyond their upper bound of applicability
        for j=1:length(indexArray)
            if (i > tzRuleVersions2[indexArray[j] + 1].toYear)
                # The rule version is no longer applicable
                if (!upperLocated)
                    #   The upper bound of the applicable rule set has been found
                    versionSet = AppRuleVersionSet(begYear, i-1, indexArray)
  #                  print(key, " D ", begYear, " ", i-1, " ", indexArray, "\n")
                    addRuleSet(ruleVersionSet,key,versionSet)
                    upperLocated = true
                    begYear = i
                end
                # Mark the rule versions for deletion
                push!(deletedIndices,j)
            end
        end
        # Delete the elements from the array in reverse order
        if (length(deletedIndices) > 0)
            for k = 1:length(deletedIndices)
                deleteat!(indexArray,deletedIndices[end+1-k])
            end
        end
        # Check for rule versions that have become applicable
        for j = firstInd:lastInd
            if (tzRuleVersions2[j].fromYear == i)
                # A new rule version has become applicable
                if (!upperLocated)
                    #   The upper bound of the applicable rule set has been found
                    #   Check whether the rule set is empty.  This can occur when there
                    #     is a gap - i.e., years with rules, followed by years with no
                    #     rules, followed by years with rules.  See Finland
                    if (length(indexArray) > 0)
                        versionSet = AppRuleVersionSet(begYear, i-1, indexArray)
  #                      print(key, " A ", begYear, " ", i-1, " ", indexArray, "\n")
                        addRuleSet(ruleVersionSet,key,versionSet)
                    end
                    begYear = i
                    upperLocated = true
                end
                # Add the rule version to the set
                # Change from 1 based array to 0 based array
                push!(indexArray,j-1)
            end
        end
    end
    # If any rules remain, they are perpetual
    if (length(indexArray) > 0)
        versionSet = AppRuleVersionSet(begYear, 1000000000, indexArray)
        addRuleSet(ruleVersionSet,key,versionSet)
#        print(key, " L ", begYear, " ", 1000000000, " ", indexArray, "\n")
    end
end
# Useful AppRuleVersionSet utilities
function getIndexArrayLength(x::AppRuleVersionSet)
    return length(x.versionIndexArray)
end
function combo(x::Array{AppRuleVersionSet})
    return maximum(map(getIndexArrayLength,x))
end
function getRuleVersionEffectiveMonth(x::UInt32)
    return tzRuleVersions2[x+1].transMonth
end
maxRuleVersions =  maximum(map(length,values(tzRules)))
numRules = length(keys(tzRules))
numRuleVersionSets = reduce(+,map(length,values(ruleVersionSet)))
maxRuleVersionSets = maximum(map(length,values(ruleVersionSet)))
maxRuleVersionsPerSet = maximum(map(combo,values(ruleVersionSet)))
print("Number of rules: ", string(numRules), "\n")
print("Maximum number of rule versions for a rule: ", string(maxRuleVersions),"\n")
print("Number of rule version sets: ", string(numRuleVersionSets),"\n")
print("Maximum number of rule version sets for a rule: ", string(maxRuleVersionSets),"\n")
print("Maximum number of rule versions in a set: ", string(maxRuleVersionsPerSet),"\n")
print("Maximum saved time (", maxTZ, "): ", string(maxSavedHours),
                              ":", string(maxSavedMinutes), "\n")
#
# Create the rules.h file
#
outfile = open(TZHeaderDir * "rules.h","w")
# Write the definition of the TZRule struct
structDef =
"#include <stdint.h>

#define maxRuleNameLength ** //  Maximum rule name length observed in IANA time zone files is %%
#define numRuleVersions xxxx // Number of rule versions in IANA time zone files
#define numRules yyyy // Number of rules in IANA time zone files
#define numRuleVersionSets zzzz // Number of rule version sets - derived from rule versions array
#define maxRuleVersionsPerSet qqqq // Maximum number of rule versions in a set
                                   //  derived from rule versions array
#define maxPeriodSegments rrrr // Maximum number of segments into which a year can be divided,
                               //  where each segment has a different s-w offset
                               //  maxRuleVersionsPerSet + 5

typedef struct
{
    // A time zone rule, created from the IANA timezone files
    char name[maxRuleNameLength]; // Rule name
    uint32_t fromYear; // First year rule is applicable.  No gigayears required as all
                       //  time zone date are from 1900 through 2100.
    uint32_t toYear;   //  Last year rule is applicable.
    uint8_t transMonth; // Month during which transiton to or from daylight savings time occurs
             						/* A rule is used to specify the day of the month during which the transitoin to
               							daylight savings time occurs.  The form of this rule is, for example,
               							the third Sunday on or after the 4th day of the month, or, generically,
               							the jth k day on or after the ith of the month, where i is afterDayOfMonth,
               							j is weekNumber and k is dayOfWeek. */
    uint8_t afterDayOfMonth; /* The day of the month the transition takes place, if weekNumber is 0,
                   								or the day on which to begin counting if weekNumber > 0.
                   								For example, if the third day of the month is a Sunday, and the
                   								rule is the first Sunday on or after the third of the month,
                   								the transition would occur on the third.  If it were the second
                   								Monday on or after the third of the month, it would occur on the
                   								11th.  */
    uint8_t weekNumber; // 1 = first week of the month, 2 = second week of the month..., 6 = last week of the month
    uint8_t dayOfWeek; // 0 = Sunday, 1 = Monday, etc.
    uint8_t hourOfTransition; // Hour of the day the transition is to take place
    uint8_t minuteOfTransition; // Minute of the day the transition is to take place
    uint8_t secondOfTransition; // Second of the day the transition is to take place
    uint8_t hundredthOfTransition; // Hundredth of a second of the day the transition is to take place
    uint8_t frameOfTransition; // 0 = universal time, 1 = local standard time, 2 = local wall time
    // Saved time is applied to the offset in the time zone applying the rule
    uint8_t hoursSaved; // Hours saved vs standard time after the transition
    uint8_t minutesSaved; // Minutes saved vs standard time after the transition
    char abbrevChars[^]; // Characters used to modify the time zone abbreviation when the rule is in effect, e.g., D or S
             						 //  maximum length observed in IANA time zone files is @
} TZRule;

typedef struct
{
  // A time zone rule set, containing the rule versions applicable for a range
  //  of years for a rule
  char ruleName[maxRuleNameLength]; // Rule name
  uint32_t fromYear; // First year the rule set is applicable
  uint32_t toYear; // Last year the rule set is applicable
  uint8_t numRuleVers; // Number of rule versions in the set
  uint32_t ruleVersionIndices[maxRuleVersionsPerSet]; // Indices of the rule versions in the set
} TZRuleSet;

static const TZRule TZRuleVersions[numRuleVersions] =
{\n"
# Replace the calculated values in structDef
structDef2 = replace(
              replace(
                replace(
                  replace(
                    replace(
                      replace(
                        replace(
                          replace(
                            replace(structDef,"**" => string(maxRNameLength + 1)),
                                "%%" => string(maxRNameLength)),
                                  "^" => string(maxAbbrevLength + 1)),
                                    "@" => string(maxAbbrevLength)),
                                      "xxxx" => string(nrules)),
                                       "yyyy" => string(numRules)),
                                        "zzzz" => string(numRuleVersionSets)),
                                         "qqqq" => string(maxRuleVersionsPerSet)),
                                          "rrrr" => string(maxRuleVersionsPerSet + 5))
# Write the beginning of the file
write(outfile, structDef2)
# Write each rule version
earliestFromYear = 1000000000
latestToYear = 0
for rule in tzRuleVersions2
  write(outfile,"    \u0022",rule.name,"\u0022, ", string(rule.fromYear), ", ",
                string(rule.toYear), ", ", string(rule.transMonth), ", ",
                string(rule.afterDayOfMonth), ", ", string(rule.weekNumber), ", ",
                string(rule.dayOfWeek), ", ", string(rule.hourOfTransition), ", ",
                string(rule.minuteOfTransition), ", ", string(rule.secondOfTransition), ", ",
                string(rule.hundredthOfTransition), ", ", string(rule.basisOfTransition), ", ",
                string(rule.hoursSaved), ", ",
                string(rule.minutesSaved), ", \u0022", rule.abbrevChars, "\u0022,",
                "\n")
  # Check for a situation where the wall time of the transition could be before
  #  the fromYear or after the toYear.  This can occur when the transition occurs on
  #  January 1st or December 31st and the basis is not wall.
  #
  global earliestFromYear = rule.fromYear < earliestFromYear ? rule.fromYear : earliestFromYear
  global latestToYear = (rule.toYear > latestToYear && rule.toYear != 1000000000) ?
   rule.toYear : latestToYear
  if (((rule.transMonth == 12 || rule.transMonth == 1) && rule.basisOfTransition !=2
      && rule.dayOfWeek != 7) || # In Dec or Jan with a dayOfWeek specification
      (rule.transMonth == 12 && rule.afterDayOfMonth == 31 && rule.basisOfTransition != 2) ||
      (rule.transMonth == 1 && rule.afterDayOfMonth == 1 && rule.basisOfTransition != 2))
      print("** Potential transition issue: ", rule.name, " Month: ", rule.transMonth,
      " Day: ", rule.afterDayOfMonth, " Basis: ", rule.basisOfTransition,
      " Day of week: ", rule.dayOfWeek, "\n")
  end
end
# Write the end of the array definition
write(outfile, "};")
print("Earliest rule 'from year': ", earliestFromYear, "\n")
print("Latest rule 'to year': ", latestToYear, "\n")
# Check for out of sequence rule versions
for key in keys(tzRules)
  val = tzRules[key]
  if (length(val) > 1)
    for i = 2:length(val)
      if (parse(Int, val[i]) != parse(Int, val[i-1]) + 1)
        print("Out of sequence rule version - Rule: ", key, " Index: ", string(i),"\n")
      end
    end
  end
end
# write the rule name array
structDef =
"\n
static const char TZRules[numRules][maxRuleNameLength] =
{\n"
write(outfile,structDef)
# Write each rule name
for key in keys(tzRules)
  write(outfile,"    \u022", key, "\u022,\n")
end
write(outfile, "};")
# Write the rule index array
structDef =
"\n
static const uint32_t TZRuleIndices[numRules][2] =
{\n"
write(outfile,structDef)
# Write each rule index
for key in keys(tzRules)
  write(outfile,"    ", tzRules[key][1], ", ", tzRules[key][end], ",\n")
end
write(outfile, "};")
# Write the rule set array
structDef =
"\n
static const TZRuleSet TZRuleVersionSets[numRuleVersionSets] =
{\n"
write(outfile,structDef)
# Write each rule set in the order of the rules array, sorted by
#   the month each rule version becomes effective
for key in keys(tzRules)
  versionSets = ruleVersionSet[key]
  for set in versionSets
    # Create the output array
    effectiveMonths = map(getRuleVersionEffectiveMonth, set.versionIndexArray)
    p = sortperm(effectiveMonths)
    # Check whether any of the versions have the same effective month
    for i = 2:length(p)
      if (p[i] == p[i-1])
        print("***** Duplicate effective month! Rule = ", key, " ", set.startYear,
          " - ", set.endYear, "\n")
      end
    end
    outArray = set.versionIndexArray[p]
    # Lengthen the index array to simplify output
    while (length(outArray) < maxRuleVersionsPerSet)
      push!(outArray, 0)
    end
    write(outfile,"    \u022", key, "\u022, ", string(set.startYear), ", ",
      string(set.endYear), ", ", string(length(set.versionIndexArray)), ", ")
    for i = 1:length(outArray)
      write(outfile, string(outArray[i]), ", ")
    end
    write(outfile, "\n")
  end
end
write(outfile, "};")
# Write the rule set index array in the same order as the rule name array
structDef =
"\n
static const uint32_t TZRuleSetIndices[numRules][2] =
{\n"
write(outfile,structDef)
firstind = 0
for key in keys(tzRules)
  lastind = firstind + length(ruleVersionSet[key]) - 1
  write(outfile,"    ",  string(firstind), ", ", string(lastind), ",\n")
  global firstind = lastind + 1
end
write(outfile, "};")
close(outfile)
# Reread each of the TZ files and create the timezone header file
# Declare the array to contain the time zones
tzArray = Array{String}(undef,3000,8)
ntz = 0
ntznorep = 0
maxNameLength = 0
maxAbbrevLength = 0
maxAbbrev = ""
#  Create a dictionary with the names for each time zone.
#      Initialie it in this section and finish in links secondOfTransition
TZNames = Dict{String,Array{String}}()
for file in TZFile
  infile = open(TZFileDir * file,"r")
  lines = readlines(infile)
  withinTZ = false
  for line in lines
    # If not already within a time zone, look for the first line of a time zone
    if (!withinTZ)
        if (length(line) > 4 && line[1:1] == "Z" && line[1:4] == "Zone")
            # This is the first line of a time zone specification
            global ntz += 1
            global ntznorep += 1
            withinTZ = true
            nsubs = 0
            # Blank elements 5 through 8, so they are not undefined
            tzArray[ntz,5:8] .= ""
            for sub in split(line)
                nsubs += 1
                # Check for maximum name length
                if (nsubs == 2)
                    global maxNameLength = max(maxNameLength, length(sub))
                    addString(TZNames,sub,sub)
                end
                # Check for maximum abberviation characters length
                if (nsubs == 5 && length(sub) > maxAbbrevLength)
                    global maxAbbrevLength = length(sub)
                    global maxAbbrev = sub
                end
                # Sometimes a comment and another word or phrase is added at the end of a time zone specification
                if (sub == "#")
                    break
                elseif (nsubs > 1)
                    tzArray[ntz, nsubs - 1] = sub
                end
            end
            # Check for 1 line time zone specification
            if (nsubs == 5)
                # No UNTIL field
                withinTZ = false
            end
        end
    else
        # Within a time zone.  Search for continuation lines
        # Every continuation line must have at least GMTOFF, RULES and FORMAT fields
        # Every time zone must end with a blank UNTIL field
        if (length(line) > 0 && line[1] != '#' && occursin(":", line))
            # It is a good assumption that a non-comment line with a colon is a time zone continuation lines
            # Consider a continuation line to be a new time zone with the same name
            ntz += 1
            nsubs = 0
            tzArray[ntz, 1] = tzArray[ntz - 1, 1]
            # Blank elements 5 through 8, so they are not undefined
            tzArray[ntz,5:8] .= ""
            for sub in split(line)
                nsubs += 1
                # Check for maximum abberviation characters length
                if (nsubs == 3 && length(sub) > maxAbbrevLength)
                    global maxAbbrevLength = length(sub)
                    global maxAbbrev = sub
                end
                # Sometimes a comment and another word or phrase is added at the end of a time zone specification
                if (sub == "#")
                    nsubs -= 1
                    break
                else
                    tzArray[ntz, nsubs + 1] = sub
                end
            end
            # Check end of time zone specification
            if (nsubs == 3)
                # No UNTIL field
                withinTZ = false
            end
        end
    end
  end
  close(infile)
end
print("Processed " * string(ntznorep) * " time zones.\n")
print("Processed " * string(ntz) * " time zone versions.\n")
print("Longest time zone name: " * string(maxNameLength) * " characters.\n")
print("Longest time zone abberviation (" * maxAbbrev * "): " * string(maxAbbrevLength) * " characters.\n")
#
# Create an array of TimeZones from the tzArray
#
timezones = Array{TimeZone}(undef, ntz)
#  Create a dictionary of time zone version indices.
TZIndices = Dict{String,Array{String}}()
# Process the time zone array
maxGMTHours = 0
maxGMTMinutes = 0
maxNegGMTHours = 0
maxNegGMTMinutes = 0
maxTZ = ""
maxNegTZ = ""
for i = 1:ntz
  # Populate the index dictionary
  addString(TZIndices,tzArray[i,1],string(i-1))
  # **** GMTOFF field ****
  # Check for a negative sign
  if (tzArray[i,2][1] == '-')
    isneg = 1
    gmtoff = tzArray[i,2][2:end]
  else
    isneg = 0
    gmtoff = tzArray[i,2]
  end
  # Default all values to 0
  hour = 0
  minute = 0
  second = 0
  hundredth = 0
  # A '-' indicates 0
  if (gmtoff[1] != '-')
    x = split(gmtoff,":")
    hour = parse(Int,x[1])
    if (length(x) > 1)
      minute = parse(Int,x[2])
      if (isneg == 0 && 60*hour + minute > 60*maxGMTHours + maxGMTMinutes)
        global maxGMTHours = hour
        global maxGMTMinutes = minute
        global maxTZ = tzArray[i,1]
      elseif (isneg == 1  && 60*hour + minute > 60*maxNegGMTHours + maxNegGMTMinutes)
        global maxNegGMTHours = hour
        global maxNegGMTMinutes = minute
        global maxNegTZ = tzArray[i,1]
      end
      if (length(x) == 3)
        # Split again to locate fractional seconds
        y = split(x[3],".")
        second = parse(Int,y[1])
        if (length(y) == 2)
          hundredth = parse(Int,y[2])
          if (hundredth < 10)
            hundredth *= 10
          end
        end
      end
    end
  end
  # Try to match the offset to one of the candidates fractional second
  #   GMT offset specifications.  Rounding isn't always done consistently
  #   in the files, so check both rounded and truncated values.
  for icand = 1:numCandidates
    if ((candidateFracSpec[icand,2] == hour && candidateFracSpec[icand,3] == minute &&
          candidateFracSpec[icand,4] == second) || (candidateFracSpec[icand,5] == hour &&
          candidateFracSpec[icand,6] == minute && candidateFracSpec[icand,7] == second))
      # Found a match
      candidateFracSpec[icand,1] += 1
      global numMatches += 1
      # Set the values to the more precise fractional second specification
      hour = candidateFracSpec[icand,5]
      minute = candidateFracSpec[icand,6]
      second = candidateFracSpec[icand,7]
      hundredth = candidateFracSpec[icand,8]
      break
    end
  end
  # **** RULES field ****
  # Default the hours and minuts saved to zero
  hourssaved = 0
  minutessaved = 0
  rule = tzArray[i,3]
  if (rule == "-")
    dstr = ""
  elseif (occursin(":",rule))
    x = split(rule,":")
    hourssaved = parse(Int,x[1])
    if (length(x) == 2)
      minutessaved = parse(Int,x[2])
    end
    dstr = ""
  else
    # Reference to a rule name
    dstr = rule
  end
  # **** FORMAT field ****
  format = (tzArray[i,4] != "-" ? tzArray[i,4] : "")
  # **** UNTIL field ****
  # Default all fields, except for year, to earliest possible value
  untilMonth = 1
  untilDayOfMonth = 1
  untilHour = 0
  untilMinute = 0
  untilSecond = 0
  untilHundredth = 0
  untilBasis = 2
  # If the field is blank, the time zone is effective until max time
  if (tzArray[i,5] == "")
    untilYear = 1000000000
  else
    # At least a year is specified
    untilYear = parse(Int,tzArray[i,5])
    # Next field, if specified, is the month
    if (tzArray[i,6] != "")
      for j in eachindex(months)
        if (tzArray[i,6] == months[j])
          untilMonth = j
          break
        end
      end
    end
    # Next field, if specified, is the day of the month
    if (tzArray[i,7] != "")
      # Sometimes, this field contains a 'lastddd' specification
      if (occursin("last",tzArray[i,7]))
        # Use date functions to resolve to an integer
        # Extract the day of the week
        for j in eachindex(daysj)
          if (tzArray[i,7][5:7] == daysj[j])
            # Find the day of the month of the last day of the until year and month
            untilDayOfMonth = Dates.dayofmonth(Dates.tolast(Dates.Date(untilYear,
                                untilMonth, 1), j))
            break
          end
        end
        # Sometimes, this field contains a 'ddd>=1' specification
      elseif (occursin(">=1",tzArray[i,7]))
        # Use date functions to resolve to an integer
        # Extract the day of the week
        for j in eachindex(daysj)
          if (tzArray[i,7][1:3] == daysj[j])
            # Find the day of the month of the last ddd of the until year and month
            untilDayOfMonth = Dates.dayofmonth(Dates.tofirst(Dates.Date(untilYear,
                                untilMonth, 1), j))
            break
          end
        end
      else
        untilDayOfMonth = parse(Int,tzArray[i,7])
      end
    end
    # Next field, if specified, contains the time and basis
    # Determin the UNTIL time basis by checking the last character of the AT field
    if (tzArray[i,8] != "")
      v = tryparse(Int,string(tzArray[i,8][end]))
      if (v == nothing)
        if (tzArray[i,8][end] == 'u')
          untilBasis = 0
        elseif (tzArray[i,8][end] == 's')
          untilBasis = 1
        else
          untilBasis = 2
        end
        # Remove the last character for subsequent processing
        untilTime = tzArray[i,8][1:end-1]
      else
        untilTime = tzArray[i,8]
      end
      # A - indicates midnight
      if (untilTime[1] != '-')
        x = split(untilTime,":")
        untilHour = parse(Int,x[1])
        if (length(x) > 1)
          untilMinute = parse(Int,x[2])
          if (length(x) == 3)
            # Split again to locate fractional seconds
            y = split(x[3],".")
            untilSecond = parse(Int,y[1])
            if (length(y) == 2)
              untilHundredth = parse(Int,y[2])
              if (untilHundredth < 10)
                untilHundredth *= 10
              end
            end
          end
        end
      end
    end
  end
  # Create a TimeZone
  timezones[i] = TimeZone(tzArray[i,1], isneg, hour, minute, second, hundredth, dstr,
                  hourssaved, minutessaved, format, untilYear, untilMonth,
                  untilDayOfMonth, untilHour, untilMinute, untilSecond,
                  untilHundredth, untilBasis)
end
function getUntilYear(x::TimeZone)
    return x.untilYear
end
print("Earliest until year: ", minimum(map(getUntilYear, timezones)), "\n")
print("Number of candidate fractional second GMTOff specifications matched: ", string(numMatches)," \n")
print("Maximum positive GMT Offset (", maxTZ, "): ", string(maxGMTHours),
                              ":", string(maxGMTMinutes), "\n")
print("Maximum negative GMT Offset = (", maxNegTZ, "): -", string(maxNegGMTHours), ":", string(maxNegGMTMinutes), "\n")
#
# Process Links
#
# Reread the files
nlinks = 0
maxLinkLength = 0
for file in TZFile
  infile = open(TZFileDir * file,"r")
  lines = readlines(infile)
  for line in lines
      # Search for lines starting with Link
      if (length(line) > 4)
          subs = split(line)
          if (subs[1] == "Link")
            # Found a link
              global nlinks += 1
              addString(TZNames,subs[2],subs[3])
              global maxLinkLength = max(maxLinkLength, length(subs[3]))
          end
      end
  end
  close(infile)
end
maxLinks =  maximum(map(length,values(TZNames))) - 1
print("Processed " * string(nlinks) * " links.\n")
print("Longest link name: " * string(maxLinkLength) * " characters.\n")
# Find the greatest number of links for a time zone
print("Most links for a time zone: " * string(maxLinks),"\n")
#
# Create the timezones.h file
#`
outfile = open(TZHeaderDir * "timezones.h","w")
# Write the definition of the TimeZone struct
structDef =
"#include <stdint.h>

#define maxTZNameLength cccc // Maximum time zone / link name length observed in IANA time zone files is qq
#define numTimeZones yyyy // Number of time zones in IANA time zone files
#define numTimeZoneVersions xxxx // Number of time zone versions in IANA time zone files
#define maxNamesPerTimeZone bbbb // Maximum number of names for a time zone observed in IANA time zone files

typedef struct
{
    // A time zone version, created from the IANA time zone files
    //  Every continuation line in the time zone file is another time zone version
    char name[maxTZNameLength]; // Time zone name
    // GMT Offset is the amount of time to add to UTC to get the standard time in this time zone
    uint8_t isNegativeGMTOffset;  // 0 = positive offset, 1 = negative offset
    uint8_t hoursOfGMTOffset;  // Number of hours of GMT offset
    uint8_t minutesOfGMTOffset;  // Number of minutes of GMT offset
    uint8_t secondsOfGMTOffset;  // Number of seconds of GMT offset
    uint8_t hundredthsOfGMTOffset;  // Number of hundredths of a second of GMT Offset
    // A time zone can either have no daylight savings time rule, 1 daylight savings time rule, or
    //  a specified time offset.
    // If the rule field is blank, use the offset.  When no rule is specified, set the offset to zero,
    //  to handle the second case, above.
    char daylightSavingsTimeRule[maxRuleNameLength];
    uint8_t hoursSaved;  // Number of hours saved vs standard time after the transition
    uint8_t minutesSaved;  // Number of minutes saved vs standard time after the transition
    // The following field contains the time zone abbreviation.  If it includes a slash,
    //  the field contains the abbreviations for both standard and daylight savings time.
    // If it includes %s, this is a place holder for the characters provided by the daylight
    //  savings time rules.
    char abbrevChars[zzz]; // Maximum length observed in IANA time zone files is @
    // The until datetime indicates the date and time until which the time zone is effective,
    //  using the GMT offset and daylight savings time rules of the time zone.
    //  The default basis of this datetime is local wall time.`
    // The basis can be modified to local standard time or to universal time.
    uint32_t untilYear;  // Year of the datetime until which the time zone is effective.  No need to consider gigayear.
    uint8_t untilMonth;  // Month of the datetime until which the time zone is effective.
    uint8_t untilDayOfMonth;  // Day of the datetime until which the time zone is effective.
    uint8_t untilHour;  // Hour of the datetime until which the time zone is effective.
    uint8_t untilMinute;  // Minute of the datetime until which the time zone is effective.
    uint8_t untilSecond;  // Second of the datetime until which the time zone is effective.
    uint8_t untilHundredth;  // Hundredth of a second of the datetime until which the time zone is effecttive.
    uint8_t untilFrame; // 0 = universal time, 1 = local standard time, 2 = local wall time
} TimeZone;

static const TimeZone TimeZoneVersions[numTimeZoneVersions] =
{\n"
# Replace the calculated values in structDef
structDef2 = replace(
              replace(
                replace(
                  replace(
                    replace(structDef,"**" => string(maxNameLength + 1)),
                              "qq" => string(max(maxNameLength, maxLinkLength))),
                                "zzz" => string(maxAbbrevLength + 1)),
                                  "@" => string(maxAbbrevLength)),
                                    "xxxx" => string(ntz))
structDef3 = replace(
               replace(
                  replace(structDef2,"cccc" => string(max(maxNameLength+1, maxLinkLength+1))),
                                      "bbbb" => string(maxLinks+1)),
                                        "yyyy" => string(ntznorep))
# Write the beginning of the file
write(outfile, structDef3)
# Write each time zone
for tz in timezones
   write(outfile,"    \u0022",tz.name,"\u0022, ", string(tz.isNegativeGMTOffset), ", ",
                string(tz.hoursOfGMTOffset), ", ", string(tz.minutesOfGMTOffset), ", ",
                string(tz.secondsOfGMTOffset), ", ", string(tz.hundredthsOfGMTOffset), ", ",
                "\u0022", tz.daylightSavingsTimeRule, "\u0022, ", string(tz.hoursSaved), ", ",
                string(tz.minutesSaved), ", \u0022", tz.abbrevChars, "\u0022, ",
                string(tz.untilYear), ", ", string(tz.untilMonth), ", ",
                string(tz.untilDayOfMonth), ", ", string(tz.untilHour), ", ",
                string(tz.untilMinute), ", ", string(tz.untilSecond), ", ",
                string(tz.untilHundredth), ", ", string(tz.untilBasis), ",\n")
    # Check for potential transition issues.  Specifically, where the Until datetime is
    #  either on 12/31 or 1/1 and the bais is not wall time.  In these cases, the until time
    #  could move to the next or previous year when converted to wall time, complicating the
    #  matching with the year of the rule set
    if (tz.untilBasis != 2 && tz.daylightSavingsTimeRule != "" &&
        ((tz.untilMonth == 12 && tz.untilDayOfMonth == 31) ||
        (tz.untilMonth == 1 && tz.untilDayOfMonth == 1)))
        print("*** Potential transition issue in time zone ", tz.name, " - month: ",
        tz.untilMonth, " day: ", tz.untilDayOfMonth, "\n")
    end
end
# Write thte end of the array
write(outfile, "};")
# Create a TZNames array in the timezone.h file
structDef =
"

static const char TimeZones[numTimeZones][maxNamesPerTimeZone][maxTZNameLength] =
{\n"
write(outfile, structDef)
# Write each link
k = collect(keys(TZNames))
for i in eachindex(k)
    v = TZNames[k[i]]
    lenThisLine = 0
    write(outfile,  "    \u0022", k[i], "\u0022, ")
    lenThisLine += 8 + length(k[i])
    for j in eachindex(v)
      if (j > 1)
        write(outfile, "\u0022", v[j], "\u0022, ")
        lenThisLine += 4 + length(v[j])
        # Limit line length
        if (lenThisLine > 100)
            write(outfile, "\n        ")
            lenThisLine = 8
        end
      end
    end
    if length(v) < maxLinks + 1
        if (lenThisLine + 4 * (maxLinks + 1 - length(v)) > 110)
            write(outfile, "\n        ")
        end
        for m = length(v) + 1 :  maxLinks + 1
            write(outfile, "\u0022\u0022, ")
        end
    end
    write(outfile, "\n")
end
# Write thte end of the array
write(outfile, "};")
structDef =
"

static const uint32_t TimeZoneIndices[numTimeZones][2] =
{\n"
write(outfile, structDef)
# Write each time zone index
for key in keys(TZIndices)
  write(outfile,"    ", TZIndices[key][1], ", ", TZIndices[key][end], ",\n")
  # Check whether a dstr has been specified for the first version of a time zone
  firstVer = parse(Int, TZIndices[key][1])
  if (timezones[firstVer+1].daylightSavingsTimeRule != "")
    print("***** First version of time zone ", key, " has a dstr.\n")
  end
end
write(outfile, "};")
close(outfile)
end # Ends module FileProcessor

max_stars_repo_path: Julia Integration/TZFileReader.jl

max_stars_repo_name: JayDweck/Ultimate-Datetime

max_issues_repo_path: Julia Integration/TZFileReader.jl

max_issues_repo_name: JayDweck/Ultimate-Datetime

max_forks_repo_path: Julia Integration/TZFileReader.jl

max_forks_repo_name: JayDweck/Ultimate-Datetime


text:
@testset "Time grid generation" begin
    @test isempty( MP.generateTimeGrid( mpSim, -12 ) )
    @test MP.generateTimeGrid( mpSim, 12 ) == [0.0]
    run( mpSim.sim )
    @test MP.generateTimeGrid( mpSim, 12 ) == collect( 0.0:12.0:300.0 )
end  # @testset "Time grid generation"

ext: jl


text:
using TopOptLevelSet
using Test

@testset "TopOptLevelSet.jl" begin
    # Write your tests here.
end


text:
using Base.BLAS

include("units.jl")

immutable Layer{T<:FloatingPoint}
	ni::Int
	no::Int
	W::Matrix{T}
	B::Vector{T}
	GW::Matrix{T}
	GB::Vector{T}
	UW::Matrix{T}
	UB::Vector{T}
	NET::Vector{T}
	ACT::Vector{T}
	DE_DNET::Vector{T}
	E::Vector{T}
	DE_DYH::Vector{T}
	DELTAS::Vector{T}
	activation::ASCIIString
	activation_function!::Function
	activation_backward!::Function

	function Layer(ni::Int, no::Int, activation::ASCIIString)
		if ni > 0 && no > 0
			# Initialize weights and biases for layer.
			sigma = 1e-2
			W = zeros(T, (ni, no))
			for i = 1:length(W)
				W[i] = sigma * randn()
			end
			B = zeros(T, no)
			for i = 1:length(B)
				B[i] = sigma * randn()
			end
			# Initialize storage for gradient information.
			GW, GB = zeros(W), zeros(B)
			UW, UB = zeros(W), zeros(B)
			# Initialize storage for upper level units.
			NET = zeros(T, no)
			ACT = zeros(T, no)
			DE_DNET = zeros(T, no)
			E = zeros(T, no)
			DE_DYH = zeros(T, no)
			DELTAS = zeros(T, ni)
			# Set activation function for layer.
			activation, activation_function!, activation_backward! = activation_function_selector(activation)
			# Create and return the object.
			new(ni, no, W, B, GW, GB, UW, UB, NET, ACT, DE_DNET, E, DE_DYH, DELTAS, activation, activation_function!, activation_backward!)
		else
			error("Invalid number of units used to initialize Layer object (ni=$ni; no=$no) to create Layer object.")
		end
	end
end

# Propagate delta quantities coming from the layer above downward through the layer.
function backward!{T<:FloatingPoint}(L::Layer{T}, DELTAS_ABOVE::Vector{T})
	@inbounds begin
		# Set DE_DNET on the layer.
		L.activation_backward!(L.NET, L.ACT, DELTAS_ABOVE, L.DE_DNET)
		# Set DELTAS on the layer.
		fill!(L.DELTAS, 0.0)
		gemv!('N', 1.0, L.W, L.DE_DNET, 1.0, L.DELTAS)
	end
end

# Propagate an incoming pattern forward through a layer.
function forward!{T<:FloatingPoint}(L::Layer{T}, IN::AbstractVector{T})
	@inbounds begin
		# Calculate net values.
		blascopy!(L.no, L.B, 1, L.NET, 1)
		gemv!('T', 1.0, L.W, IN, 1.0, L.NET)
		# Set activations and gradient information related to activations.
		L.activation_function!(L.NET, L.ACT)
	end
end


text:
"""
Abstract type for Service Formulations (a.k.a Models)

# Example
```julia
import PowerSimulations
const PSI = PowerSimulations
struct MyServiceFormulation <: PSI.AbstractServiceFormulation
```
"""
abstract type AbstractServiceFormulation end
abstract type AbstractReservesFormulation <: AbstractServiceFormulation end

function _check_service_formulation(
    ::Type{D},
) where {D <: Union{AbstractServiceFormulation, PSY.Service}}
    if !isconcretetype(D)
        throw(
            ArgumentError(
                "The service model must contain only concrete types, $(D) is an Abstract Type",
            ),
        )
    end
end

"""
Establishes the model for a particular services specified by type. Uses the keyword argument
`use_service_name` to assign the model to a service with the same name as the label in the
template. Uses the keyword argument feedforward to enable passing values between operation
model at simulation time

# Arguments
-`::Type{D}`: Power System Service Type
-`::Type{B}`: Abstract Service Formulation

# Accepted Key Words
- `feedforward::Array{<:AbstractAffectFeedForward}` : use to pass parameters between models
- `use_service_name::Bool` : use the label as the name for the service

# Example
```julia
reserves = ServiceModel(PSY.VariableReserve{PSY.ReserveUp}, RangeReserve)
```
"""
mutable struct ServiceModel{D <: PSY.Service, B <: AbstractServiceFormulation}
    component_type::Type{D}
    formulation::Type{B}
    feedforward::Union{Nothing, AbstractAffectFeedForward}
    use_service_name::Bool
    function ServiceModel(
        ::Type{D},
        ::Type{B};
        feedforward::Union{Nothing, AbstractAffectFeedForward} = nothing,
        use_service_name::Bool = false,
    ) where {D <: PSY.Service, B <: AbstractServiceFormulation}
        _check_service_formulation(D)
        _check_service_formulation(B)
        new{D, B}(D, B, feedforward, use_service_name)
    end
end

get_component_type(m::ServiceModel) = m.component_type
get_formulation(m::ServiceModel) = m.formulation
get_feedforward(m::ServiceModel) = m.feedforward

function _set_model!(dict::Dict, key::Tuple{String, Symbol}, model::ServiceModel)
    if haskey(dict, key)
        @info("Overwriting $(key) existing model")
    end
    dict[key] = model
    return
end

function _set_model!(
    dict::Dict,
    service_name::String,
    model::ServiceModel{D, B},
) where {D <: PSY.Service, B <: AbstractServiceFormulation}
    if !model.use_service_name
        throw(
            IS.ConflictingInputsError(
                "The model provided has use_service_name false. This method can't be used",
            ),
        )
    end
    _set_model!(dict, (service_name, Symbol(D)), model)
    return
end

function _set_model!(
    dict::Dict,
    model::ServiceModel{D, B},
) where {D <: PSY.Service, B <: AbstractServiceFormulation}
    if model.use_service_name
        throw(
            IS.ConflictingInputsError(
                "The model provided has use_service_name set to true and no service name was provided. This method can't be used",
            ),
        )
    end
    _set_model!(dict, (NO_SERVICE_NAME_PROVIDED, Symbol(D)), model)
    return
end


text:
input = 324

# Part 1
a = [0]
i = 0
for t in 1:2017
    i = (i + input) % t + 1
    insert!(a, i + 1, t)
end
println(a[i + 2])

# Part 2
i = 0
result = 0
for t in 1:50_000_000
    i = (i + input) % t + 1
    if i == 1
        result = t
    end
end
println(result)


text:
using NiLang
using Test

@testset "identity" begin
    x, y = 0.2, 0.5
    @instr x += identity(y)
    @test x == 0.7 && y==0.5
end

@testset "*, /" begin
    x, y, out = 2.0, 2.0, 1.0
    @instr out += x * y
    @test x == 2.0 && y == 2.0 && out == 5.0
    x, y, out = 2.0, 2.0, 1.0
    @instr out += x / y
    @test x == 2.0 && y == 2.0 && out == 2.0

    out = Fixed43(0.0)
    x = 1
    @instr out += x/2
    @test out === Fixed43(0.5)
    @instr out -= x/2
    @test out === Fixed43(0.0)
end

@testset "SWAP" begin
    x, y = 1, 2
    @instr SWAP(x, y)
    @test x == 2 && y == 1
end

@testset "-" begin
    x = 0.3
    @instr -(x)
    @test x == -0.3
    @test check_inv(-, (x,))
end

@testset "FLIP" begin
    x = false
    @instr FLIP(x)
    @test x == true
    @test check_inv(FLIP, (x,))
end

@testset "ROT" begin
    x, y, θ = 0.0, 1.0, π
    @test check_inv(ROT, (x, y, θ); verbose=true)
    @test check_inv(IROT, (x, y, θ); verbose=true)
end

@testset "INC, DEC" begin
    x = Int32(2)
    @instr INC(x)
    @test x === Int32(3)
    @instr DEC(x)
    @test x === Int32(2)
    @instr (x,) |> INC |> DEC |> DEC
    @test x === Int32(1)
end


text:
module TestDynamicExamples

import Modia3D
using  Modia3D.Test
import Modia3D.ModiaMath

# Test Modia3D - Dynamic Examples

# examples/dynamics
dynamicsPath = joinpath(Modia3D.path, "examples", "dynamics")

include(joinpath(dynamicsPath, "Simulate_Pendulum.jl"))
import .Simulate_Pendulum

include(joinpath(dynamicsPath, "Simulate_DoublePendulum.jl"))
import .Simulate_DoublePendulum

@testset "Modia3D: Simulate_Pendulum.jl" begin
   pendulum = Modia3D.SimulationModel( Simulate_Pendulum.Pendulum(Lx=1.6, m=0.5) )
   result = ModiaMath.simulate!(pendulum, stopTime=5.0, tolerance=1e-8, log=false)
   phi = result.series["rev.phi"]
   w   = result.series["rev.w"]

   @test isapprox(phi[end], -2.61422; atol=1e-5 )
   @test isapprox(w[end]  , -3.51319; atol=1e-5 )
end

@testset "Modia3D: Simulate_DoublePendulum.jl" begin
   doublePendulum = Modia3D.SimulationModel( Simulate_DoublePendulum.DoublePendulum(Lx=1.0, m=1.0) )
   result = ModiaMath.simulate!(doublePendulum, stopTime=5.0, tolerance=1e-8, log=false)
   phi1 = result.series["rev1.phi"]
   phi2 = result.series["rev2.phi"]
   w1   = result.series["rev1.w"]
   w2   = result.series["rev2.w"]

   @test isapprox(phi1[end], -0.211087; atol=1e-3 )
   @test isapprox(phi2[end],  0.458381; atol=1e-3 )
   @test isapprox(w1[end]  ,  1.51728 ; atol=1e-2 )
   @test isapprox(w2[end]  , -0.928819; atol=1e-2 )
end

end


text:
jlpath = readchomp(`which julia`)

mkpath("../julia")
run(`ln -s $jlpath ../julia/julia`)

Pkg.add("Atom")

if get(ENV, "ATOMJL", "") == "master"
  Pkg.checkout("Atom")
end

using Atom


text:
include(joinpath(@__DIR__, "../../src/FLS.jl"))

"""
Permutations of a set, i.e. invertible endo-functions.
"""
permschema = @acset LabeledGraph begin
    V = 1
    E = 2
    vlabel = [:X]
    elabel = [:f, :f⁻¹]
    src = [1,1]
    tgt = [1,1]
end

perm = FLS(:perm, permschema, [], [(:inv, [:f, :f⁻¹], [])])


text:
using HDF5

Base.Base.@kwdef mutable struct BeamInfo <: Info
    ras::AbstractVector{Float64} = Vector{Float64}(undef, 0)
    decs::AbstractVector{Float64} = Vector{Float64}(undef, 0)
    src_names::AbstractVector{String} = Vector{String}(undef, 0)
end # mutable struct BeamInfo

#=
function to_hdf5(group::HDF5.Group, info::BeamInfo)
    invoke(to_hdf5, Tuple{HDF5.Group, Info}, group, info)
end
=#


text:
abstract type AbstractOctNode end

"""
Struct to hold an OctTree node

* `subnodes`: matrix of subnodes
* `father`: `AbstractOctNode` or `Nothing`, for the node from which the object branches
* `lims`: geometric limits of the node in each axis (shape `(N, 2)`)
* `adjacent`: matrix of `AbstractOctNode` objects referring to `OctTree` nodes and leaves,
    with the same shape as `lims`
"""
mutable struct OctNode{N} <: AbstractOctNode
    subnodes::Array{AbstractOctNode, N}
    father::Union{AbstractOctNode, Nothing}
    lims::Matrix{Float64}
    adjacent::Matrix{Union{AbstractOctNode, Nothing}}
end

"""
Struct to hold an OctTree leaf

* `ind`: index for a leaf's corresponding degree of freedom. Defaults to zero in constructors,
    and is set by `set_indexing!`.
    Can also refer to an array of leaf indices resulting from the merging of a branch
* `father`: `AbstractOctNode` or `Nothing`, for the node from which the leaf branches
* `lims`: geometric limits of the node in each axis (shape `(N, 2)`):
* `adjacent`: matrix of `AbstractOctNode` objects referring to `OctTree` nodes and leaves,
    with the same shape as `lims`
"""
mutable struct OctLeaf <: AbstractOctNode
    ind::Union{Int64, AbstractArray}
    father::Union{AbstractOctNode, Nothing}
    lims::Matrix{Float64}
    adjacent::Matrix{Union{AbstractOctNode, Nothing}}
end

"""
Constructor for OctLeaf

* `lims`: geometric limits for leaf in each axis(shape `(ndims, 2)`)
"""
OctLeaf(lims::Matrix{Float64}) = OctLeaf(
    0,
    nothing,
    lims,
    Union{AbstractOctNode, Nothing}[
        nothing for i=1:size(lims, 1), j=1:2
    ]
)

"""
Get area for leaf

* `lf`: leaf at hand
* `dim`: dimension to get area in
"""
function leaf_area(lf::AbstractOctNode, dim::Int64)

    a=1.0

    for i=1:size(lf.lims, 1)
        if i!=dim
            a*=(lf.lims[i, 2]-lf.lims[i, 1])
        end
    end

    a

end

"""
Get volume for leaf

* `lf`: leaf at hand
"""
function leaf_volume(lf::AbstractOctNode)

    v=1.0

    for i=1:size(lf.lims, 1)
        v*=(lf.lims[i, 2]-lf.lims[i, 1])
    end

    v

end

"""
Check if a given position is part of a node

* `pt`: vector of floats with position
* `node`: `AbstractOctNode`
* `level`: integer, `>=0`. Indicates how many more levels should be fetched

* return: subdomain if succeeded, nothing if failed
"""
function get_subdomain(pt::Vector{Float64}, node::AbstractOctNode, level::Int64)

    if isa(node, OctLeaf)
        for ndim=1:size(node.lims, 1)
            if pt[ndim]<node.lims[ndim, 1] || pt[ndim]>node.lims[ndim, 2]
                return nothing
            end
        end

        return node
    end

    if isa(node, OctNode)
        if level==0
            for ndim=1:size(node.lims, 1)
                if pt[ndim]<node.lims[ndim, 1] || pt[ndim]>node.lims[ndim, 2]
                    return nothing
                end
            end

            return node
        else
            for sd in node.subnodes
                ret=get_subdomain(pt, sd, level-1)

                if !isnothing(ret)
                    return ret
                end
            end

            fth=node.father

            if isnothing(fth)
                return nothing
            else
                return get_subdomain(pt, fth, level+1)
            end
        end
    end

end

"""
Obtain subdomain within which a point is inserted, using a recursive search reaching
any level of discretization

* `pt`: vector of floats with position
* `node`: `AbstractOctNode` from which to search

* return: `nothing` if not found, `AbstractOctNode` for the subdomain if found
"""
get_subdomain(pt::Vector{Float64}, node::AbstractOctNode) = get_subdomain(pt, node, typemax(Int64))

"""
Set adjacent attribute in `AbstractOctNode`
"""
function set_adjacent!(node::AbstractOctNode; recursive::Bool=false)

    fth=node.father

    if isnothing(fth)
        node.adjacent.=nothing

        return
    end

    ndim=size(node.lims, 1)
    coords=zeros(Float64, ndim)

    for cinds in CartesianIndices(node.adjacent)
        nd=cinds[1]
        naux=cinds[2]

        micro_offset=(naux==1 ? -1e-10 : 1e-10)

        for c=1:ndim
            if c!=nd
                coords[c]=(node.lims[c, 1]+node.lims[c, 2])/2
            else
                coords[c]=node.lims[c, naux]+micro_offset
            end
        end

        node.adjacent[cinds]=get_subdomain(coords, fth, 1)
    end

    if recursive && isa(node, OctNode)
        for sd in node.subnodes
            set_adjacent!(sd; recursive=true)
        end
    end

end

"""
Fragment a leaf into an `OctNode` and return the result

* `lf`: `OctLeaf` struct
"""
function fracture(lf::OctLeaf)

    nd=size(lf.lims, 1)

    original_halfdims=[(lf.lims[i, 2]-lf.lims[i, 1])/2 for i=1:nd]

    subdoms=Array{AbstractOctNode, nd}(undef, [2 for i=1:nd]...)

    newdom=OctNode{nd}(
        subdoms,
        lf.father,
        lf.lims,
        lf.adjacent
    )

    for cind in CartesianIndices(subdoms)
        lims=copy(lf.lims)

        for i=1:nd
            lims[i, 1]=(cind[i]==1 ? lf.lims[i, 1] : lf.lims[i, 1]+original_halfdims[i])
            lims[i, 2]=(cind[i]==2 ? lf.lims[i, 2] : lf.lims[i, 1]+original_halfdims[i])
        end

        newdom.subnodes[cind]=OctLeaf(
            lf.ind,
            newdom,
            lims,
            Matrix{Union{AbstractOctNode, Nothing}}(undef, nd, 2),
        )
    end

    for sd in newdom.subnodes
        set_adjacent!(sd)
    end

    newdom

end

"""
Function to fracture a cell and re-adjust the connectivity map

* `oct`: `AbstractOctNode` for an `OctNode`, father to the desired `OctLeaf`,
    indexed by `inds`. If `inds` is absent, `oct` is expected to be a leaf itself,
    and no update of the nearby cells' adjacency maps is done
* `inds`: indices to a leaf in `oct`. If absent, `oct` is assumed to be a leaf itself, but the
    adjacency maps aren't updated for nearby cells

* return: a new `OctNode` to replace `oct`
"""
function fracture!(oct::AbstractOctNode; inds::Union{Nothing, Vector{Int64}}=nothing)

    if isnothing(inds) && !isa(oct, OctLeaf)
        throw(
            error(
                "fracture!:index argument not provided, but octnode requested for fracture is not a leaf"
            )
        )
    end

    if !isnothing(inds) && isa(oct, OctLeaf)
        throw(
            error(
                "fracture!:index argument provided, but oct is a leaf"
            )
        )
    end

    if isnothing(inds)
        fr=fracture(oct)

        for adj in fr.adjacent
            if !isnothing(adj)
                set_adjacent!(adj; recursive=true)
            end
        end

        return fr
    end

    oct.subnodes[inds...]=fracture(oct.subnodes[inds...])

    for adj in oct.subnodes[inds...].adjacent
        if !isnothing(adj)
            set_adjacent!(adj; recursive=true)
        end
    end

    oct

end

"""
Merge the leaves in a branch

* `nd`: `OctNode` struct for the branch

* return: an `OctLeaf` struct for the resulting branch
"""
function merge_branch(branch::OctNode{N}) where N

    inds=[
        begin
            @assert isa(sd, OctLeaf) "merge_branch:ERROR:attempting to merge a branch with more than one downstream ramification"

            sd.ind
        end for sd in branch.subnodes
    ]

    lf=OctLeaf(
        inds,
        branch.father,
        branch.lims,
        branch.adjacent
    )

    lf

end

"""
Function to merge a branch and re-adjust the connectivity map

* `oct`: `AbstractOctNode` for an `OctNode`, father to the desired branch,
    indexed by `inds`. If `inds` is absent, `oct` is expected to be the branch itself,
    and no update of the nearby cells' adjacency maps is done
* `inds`: indices to a branch in `oct`. If absent, `oct` is assumed to be the branch itself, but the
    adjacency maps aren't updated for nearby cells

* return: a new `OctNode` to replace `oct`
"""
function merge_branch!(oct::OctNode{N}; inds::Union{Nothing, AbstractVector}=nothing) where N

    if isa(inds, AbstractVector)
        @assert N==length(inds) "merge_branch!:DEBUG:wrong dimensionality passed for indexing vector"
    end

    if isnothing(inds)
        return merge_branch(oct)
    end

    oct.subnodes[inds...]=merge_branch(oct.subnodes[inds...])

    for adj in oct.subnodes[inds...].adjacent
        if !isnothing(adj)
            set_adjacent!(adj; recursive=true)
        end
    end

    oct

end

function _rec_nleaves(t::AbstractOctNode, alrdy::Int64)

    if isa(t, OctLeaf)
        return alrdy+1
    end

    n=alrdy

    for s in t.subnodes
        n=_rec_nleaves(s, n)
    end

    n

end

"""
Get number of leaves in octree (equivalent to the number of degrees of freedom for each
variable in the solved differential equation)
"""
nleaves(t::AbstractOctNode) = _rec_nleaves(t, 0)

function _rec_setinds!(n, oldvs, newvs, upto)

    if isa(n, OctLeaf)
        upto+=1

        if isa(n.ind, AbstractArray)
            for (o, nw) in zip(oldvs, newvs)
                nw[upto]=sum(o[n.ind])/length(n.ind)
            end
        else
            for (o, nw) in zip(oldvs, newvs)
                nw[upto]=o[n.ind]
            end
        end

        n.ind=upto
    else
        for sd in n.subnodes
            upto=_rec_setinds!(sd, oldvs, newvs, upto)
        end
    end

    upto

end

"""
Function to set indices at leaves of the octree

* `msh`: `AbstractOctNode` corresponding to the first layer of resolution
* `oldvs`: arrays with values for variables as previously stocked
"""
function set_indices!(msh::AbstractOctNode, oldvs::AbstractVector...)

    nlvs=nleaves(msh)

    newvs=[zeros(Float64, nlvs) for i=1:length(oldvs)]

    ret=_rec_setinds!(msh, collect(oldvs), newvs, 0)

    @assert ret==nlvs

    newvs

end

export OctLeaf, OctNode, fracture!, merge_branch!, nleaves, set_indices!, get_subdomain


text:
using Allnix: axpy
using Allnix.Test: init!
using BenchmarkTools
using Base.Test

function sample_count()
    return Int64(100)
end

function run_axpy{T}(n::Int64, α::T)
#    t = typeof(α)
    println("Testing axpy with typeof $(typeof(α))")
    x = Vector{T}(n)
    y = Vector{T}(n)

    # - Calculate the answer - #
    # - I want to double check the calculation is actually
    # - done in the benchmark code.
    init!(x, y)
    ans_1 = α*x[1] + y[1]
    ans_n = α*x[n] + y[n]

    b = @benchmarkable axpy($α, $x, $y) setup=(init!($x,$y))
    t = BenchmarkTools.run(b, samples = sample_count())
    show(STDOUT, "text/plain", t)
    print("\n\n")

    @test ans_1 == y[1]
    @test ans_n == y[n]
end

function axpy_ifelse(α::Float64, x::Vector{Float64}, y::Vector{Float64})
    @inbounds for i in 1:length(x)
        if x[i] < 0.5
            y[i] = α*x[i] + y[i]
        else
            y[i] = α*x[i] - y[i]
        end
    end
end

function run_axpy_ifelse(n::Int64, α::Float64)
    println("Testing axpy with if-else")
    x = Vector{Float64}(n)
    y = Vector{Float64}(n)

    # - Calculate the answer - #
    init!(x, y)
    if x[1] < 0.5
        ans_1 = α*x[1] + y[1]
    else
        ans_1 = α*x[1] - y[1]
    end
    if x[n] < 0.5
        ans_n = α*x[n] + y[n]
    else
        ans_n = α*x[n] - y[n]
    end

    b = @benchmarkable axpy_ifelse($α, $x, $y) setup(init!($x, $y))
    t = BenchmarkTools.run(b, samples = sample_count())
    show(STDOUT, "text/plain", t)
    print("\n\n")

    @test ans_1 == y[1]
    @test ans_n == y[n]
end

function run()
    n::Int64 = 1000000
    run_axpy(n, 0.5)
    run_axpy(n, Float32(0.5))
    run_axpy(n, 2)
    run_axpy(n, Int32(2))

    run_axpy_ifelse(n, 0.5)
end

run()


text:
for p in ("ArgParse","JLD","Knet")
    Pkg.installed(p) == nothing && Pkg.add(p)
end

"""

This example implements an LSTM network for training and testing
character-level language models inspired by ["The Unreasonable
Effectiveness of Recurrent Neural
Networks"](http://karpathy.github.io/2015/05/21/rnn-effectiveness)
from Andrej Karpathy's blog.  The model can be trained with different
genres of text, and can be used to generate original text in the same
style.

Example usage:

* `julia charlm.jl`: trains a model using its own code.

* `julia charlm.jl --data foo.txt`: uses foo.txt to train instead.

* `julia charlm.jl --data foo.txt bar.txt`: uses foo.txt for training
  and bar.txt for validation.  Any number of files can be specified,
  the first two will be used for training and validation, the rest for
  testing.

* `julia charlm.jl --best foo.jld --save bar.jld`: saves the best
  model (according to validation set) to foo.jld, last model to
  bar.jld.

* `julia charlm.jl --load foo.jld --generate 1000 --sresult generated.txt`:
  generates 1000 characters from the model in foo.jld and saves it to
  generated.txt.

* `julia charlm.jl --help`: describes all available options.

"""
module CharLM
using Knet,AutoGrad,ArgParse,JLD

# LSTM implementation with a single matrix multiplication with
# instances in rows rather than columns.  Julia is column major, so
# horizontal concatenation and column based slicing are contiguous and
# more efficient compared to vertical concatenation and row
# slicing. In this implementation I wanted to perform a single matrix
# multiplication for all four gates rather than four (or eight)
# separate matrix multiplications for performance. Thus I concatenate
# the input and the hidden, then slice out the four gates.  Both
# operations are more efficient if instances are in rows rather than
# columns.

function lstm(weight,bias,hidden,cell,input)
    gates   = hcat(input,hidden) * weight .+ bias
    hsize   = size(hidden,2)
    forget  = sigm.(gates[:,1:hsize])
    ingate  = sigm.(gates[:,1+hsize:2hsize])
    outgate = sigm.(gates[:,1+2hsize:3hsize])
    change  = tanh.(gates[:,1+3hsize:end])
    cell    = cell .* forget + ingate .* change
    hidden  = outgate .* tanh.(cell)
    return (hidden,cell)
end

# generate model parameters for k=1:length(hidden) lstm layers
# instances are in rows, vectors are row vectors
# model[2k-1]: weight matrix for the k'th lstm layer
# model[2k]: bias vector for the k'th lstm layer
# model[end-2]: embedding matrix
# model[end-1,end]: weight and bias for final prediction
function initmodel(atype, hidden, vocab, embed)
    init(d...)=atype(xavier(d...))
    bias(d...)=atype(zeros(d...))
    model = Array{Any}(2*length(hidden)+3)
    X = embed
    for k = 1:length(hidden)
        H = hidden[k]
        model[2k-1] = init(X+H, 4H)
        model[2k] = bias(1, 4H)
        model[2k][1:H] = 1 # forget gate bias = 1
        X = H
    end
    model[end-2] = init(vocab,embed)
    model[end-1] = init(hidden[end],vocab)
    model[end] = bias(1,vocab)
    return model
end

# state[2k-1]: hidden for the k'th lstm layer
# state[2k]: cell for the k'th lstm layer
let blank = nothing; global initstate
function initstate(model, batch)
    nlayers = div(length(model)-3,2)
    state = Array{Any}(2*nlayers)
    for k = 1:nlayers
        bias = model[2k]
        hidden = div(length(bias),4)
        if typeof(blank)!=typeof(bias) || size(blank)!=(batch,hidden)
            blank = fill!(similar(bias, batch, hidden),0)
        end
        state[2k-1] = state[2k] = blank
    end
    return state
end
end

# initoptim creates optimization parameters for each numeric weight
# array in the model.  This should work for a model consisting of any
# combination of tuple/array/dict.
initoptim{T<:Number}(::KnetArray{T},otype)=eval(parse(otype))
initoptim{T<:Number}(::Array{T},otype)=eval(parse(otype))
# TODO: This breaks Julia4 parser:
# initoptim(a::Associative,otype)=Dict(k=>initoptim(v,otype) for (k,v) in a)
initoptim(a,otype)=map(x->initoptim(x,otype), a)

# input: Dense token-minibatch input
# returns (B,H) hidden output and newstate
function predict(model, state, input; pdrop=0)
    nlayers = div(length(model)-3,2)
    newstate = similar(state)
    for k = 1:nlayers
        input = dropout(input, pdrop)
        (newstate[2k-1],newstate[2k]) = lstm(model[2k-1],model[2k],state[2k-1],state[2k],input)
        input = newstate[2k-1]
    end
    return input,newstate
end

function generate(model, tok2int, nchar)
    int2tok = Array{Char}(length(tok2int))
    for (k,v) in tok2int; int2tok[v] = k; end
    input = tok2int[' ']
    state = initstate(model, 1)
    # Open file for saving
    if o[:sresult] != nothing
        f = open(o[:sresult],"w")
    end
    for t in 1:nchar
        embed = model[end-2][[input],:]
        ypred,state = predict(model,state,embed)
        ypred = ypred * model[end-1] .+ model[end]
        input = sample(exp.(logp(ypred)))
        print(int2tok[input])
        # Save character to file
        if o[:sresult] != nothing
            write(f, int2tok[input])
        end
    end
    # Close file if opened
    if o[:sresult] != nothing
        close(f)
    end
    println()
end

# sequence[t]: Vector{Int} token-minibatch input at time t
function loss(model, state, sequence, range=1:length(sequence)-1; newstate=nothing, pdrop=0)
    preds = []
    for t in range
        input = model[end-2][sequence[t],:]
        pred,state = predict(model,state,input; pdrop=pdrop)
        push!(preds,pred)
    end
    if newstate != nothing
        copy!(newstate, map(AutoGrad.getval,state))
    end
    pred0 = vcat(preds...)
    pred1 = dropout(pred0,pdrop)
    pred2 = pred1 * model[end-1]
    pred3 = pred2 .+ model[end]
    logp1 = logp(pred3,2)
    nrows,ncols = size(pred3)
    golds = vcat(sequence[range[1]+1:range[end]+1]...)
    index = similar(golds)
    @inbounds for i=1:length(golds)
        index[i] = i + (golds[i]-1)*nrows
    end
    # pred3 = Array(pred3) #TODO: FIX BUGGY REDUCTION CODE FOR KNETARRAY IF PRED3 TOO BIG
    logp2 = logp1[index]
    logp3 = sum(logp2)
    return -logp3 / length(golds)
end

lossgradient = grad(loss)

function avgloss(model, sequence, S)
    T = length(sequence)
    B = length(sequence[1])
    state = initstate(model, B)
    total = count = 0
    for i in 1:S:T-1
        j = min(i+S-1,T-1)
        n = j-i+1
        total += n * loss(model, state, sequence, i:j; newstate=state)
        count += n
    end
    return total / count
end

function bptt(model, sequence, optim, S; pdrop=0)
    T = length(sequence)
    B = length(sequence[1])
    state = initstate(model, B)
    for i in 1:S:T-1
        j = min(i+S-1,T-1)
        grads = lossgradient(model, state, sequence, i:j; newstate=state, pdrop=pdrop)
        update!(model, grads, optim)
    end
end

# It turns out convergence is much faster if we keep bptt length small
# initially. So we take user's o[:seqlength] as the upper limit and
# bound the bptt length by the epoch number.
function train!(model, data, tok2int, o)
    global optim = initoptim(model,o[:optimization])
    if o[:fast]
        for epoch=1:o[:epochs]
            bptt(model, data[1], optim, min(epoch,o[:seqlength]); pdrop=o[:dropout])
        end
        return
    end
    report(ep)=(l=map(d->avgloss(model,d,100), data);println((:epoch,ep,:loss,l...));l)
    @time losses = report(0)
    devset = ifelse(length(data) > 1, 2, 1)
    devlast = devbest = losses[devset]
    for epoch=1:o[:epochs]
        @time bptt(model, data[1], optim, min(epoch,o[:seqlength]); pdrop=o[:dropout])
        @time losses = report(epoch)
        if o[:gcheck] > 0
            gradcheck(loss, model, data[1], 1:min(o[:seqlength],length(data[1])-1); gcheck=o[:gcheck], verbose=true)
        end
        devloss = losses[devset]
        if devloss < devbest
            devbest = devloss
            if o[:bestfile] != nothing
                info("Saving best model to $(o[:bestfile])")
                save(o[:bestfile], "model", model, "vocab", tok2int)
            end
        end
        devlast = devloss
    end
    if o[:savefile] != nothing
        info("Saving final model to $(o[:savefile])")
        save(o[:savefile], "model", model, "vocab", tok2int)
    end
end

function minibatch(chars, tok2int, batch_size)
    chars = collect(chars)
    nbatch = div(length(chars), batch_size)
    data = [ zeros(Int,batch_size) for i=1:nbatch ]
    for n = 1:nbatch
        for b = 1:batch_size
            char = chars[(b-1)*nbatch + n]
            data[n][b] = tok2int[char]
        end
    end
    return data
end

function sample(p)
    p = convert(Array,p)
    r = rand()
    for c = 1:length(p)
        r -= p[c]
        r < 0 && return c
    end
end

function main(args=ARGS)
    global model, text, data, tok2int, o
    s = ArgParseSettings()
    s.description="charlm.jl: Character level language model based on http://karpathy.github.io/2015/05/21/rnn-effectiveness. (c) Emre Yolcu, Deniz Yuret, 2017."
    s.exc_handler=ArgParse.debug_handler
    @add_arg_table s begin
        ("--datafiles"; nargs='+'; help="If provided, use first file for training, second for dev, others for test.")
        ("--loadfile"; help="Initialize model from file")
        ("--savefile"; help="Save final model to file")
        ("--bestfile"; help="Save best model to file")
        ("--generate"; arg_type=Int; default=0; help="If non-zero generate given number of characters.")
        ("--epochs"; arg_type=Int; default=1; help="Number of epochs for training.")
        ("--hidden"; nargs='+'; arg_type=Int; default=[334]; help="Sizes of one or more LSTM layers.")
        ("--embed"; arg_type=Int; default=168; help="Size of the embedding vector.")
        ("--batchsize"; arg_type=Int; default=256; help="Number of sequences to train on in parallel.")
        ("--seqlength"; arg_type=Int; default=100; help="Maximum number of steps to unroll the network for bptt. Initial epochs will use the epoch number as bptt length for faster convergence.")
        ("--optimization"; default="Adam()"; help="Optimization algorithm and parameters.")
        ("--gcheck"; arg_type=Int; default=0; help="Check N random gradients.")
        ("--seed"; arg_type=Int; default=-1; help="Random number seed.")
        ("--atype"; default=(gpu()>=0 ? "KnetArray{Float32}" : "Array{Float32}"); help="array type: Array for cpu, KnetArray for gpu")
        ("--fast"; action=:store_true; help="skip loss printing for faster run")
        ("--dropout"; arg_type=Float64; default=0.0; help="Dropout probability.")
        ("--sresult"; help = "Save generated text to file" )
    end
    isa(args, AbstractString) && (args=split(args))
    if in("--help", args) || in("-h", args)
        ArgParse.show_help(s; exit_when_done=false)
        return
    end
    o = parse_args(args, s; as_symbols=true)
    if !o[:fast]
        println(s.description)
        println("opts=",[(k,v) for (k,v) in o]...)
    end
    o[:seed] > 0 && srand(o[:seed])
    o[:atype] = eval(parse(o[:atype]))

    # we initialize a model from loadfile, train using datafiles (both optional).
    # if the user specifies neither, train a model using the charlm.jl source code.
    isempty(o[:datafiles]) && o[:loadfile]==nothing && push!(o[:datafiles],@__FILE__) # shakespeare()

    # read text and report lengths
    text = map((@compat readstring), o[:datafiles])
    !isempty(text) && !o[:fast] && info("Chars read: $(map((f,c)->(basename(f),length(c)),o[:datafiles],text))")

    # tok2int (char_to_index) comes from the initial model if there is one, otherwise from the datafiles.
    # if there is an initial model make sure the data has no new vocab
    if o[:loadfile]==nothing
        tok2int = Dict{Char,Int}()
        for t in text, c in t; get!(tok2int, c, 1+length(tok2int)); end
        model = initmodel(o[:atype], o[:hidden], length(tok2int), o[:embed])
    else
        info("Loading model from $(o[:loadfile])")
        tok2int = load(o[:loadfile], "vocab")
        for t in text, c in t; haskey(tok2int, c) || error("Unknown char $c"); end
        model = map(p->convert(o[:atype],p), load(o[:loadfile], "model"))
    end
    !o[:fast] && info("$(length(tok2int)) unique chars.")
    if !isempty(text)
        data = map(t->minibatch(t, tok2int, o[:batchsize]), text)
        train!(model, data, tok2int, o)
    end
    if o[:generate] > 0
        generate(model, tok2int, o[:generate])
    end
    return model
end

# This allows both non-interactive (shell command) and interactive calls like:
# $ julia charlm.jl --epochs 10
# julia> CharLM.main("--epochs 10")
PROGRAM_FILE=="charlm.jl" && main(ARGS)

end  # module

# Note: 10.txt used in the sample runs below was generated using
#   head -10000 100.txt > 10.txt
# where 100.txt is the file downloaded by shakespeare().

### SAMPLE RUN 74a2e6c+ Mon Sep 19 14:03:10 EEST 2016
### Implemented multi-layer.  Removed the keepstate option fixing it to true.
### Note that winit default changed so I specify it below for comparison.
### The slight difference is due to keepstate.

# julia> CharLM.main("--data 10.txt --winit 0.3 --fast")
# charlm.jl (c) Emre Yolcu, Deniz Yuret, 2016. Character level language model based on http://karpathy.github.io/2015/05/21/rnn-effectiveness.
# opts=(:lr,4.0)(:atype,"KnetArray{Float32}")(:winit,0.3)(:savefile,nothing)(:loadfile,nothing)(:generate,0)(:bestfile,nothing)(:gclip,3.0)(:hidden,[256])(:epochs,3)(:decay,0.9)(:gcheck,0)(:seqlength,100)(:seed,42)(:embed,256)(:batchsize,128)(:datafiles,Any["10.txt"])(:fast,true)
# INFO: Chars read: [("10.txt",425808)]
# INFO: 87 unique chars.
#   1.406687 seconds (1.74 M allocations: 196.799 MB, 2.49% gc time)
# (:epoch,0,:loss,6.1075509900258)
#   4.002638 seconds (6.12 M allocations: 374.188 MB, 2.41% gc time)
#   3.990772 seconds (6.11 M allocations: 374.129 MB, 2.44% gc time)
#   4.006249 seconds (6.12 M allocations: 374.240 MB, 2.53% gc time)
#   1.405878 seconds (1.75 M allocations: 197.059 MB, 2.56% gc time)
# (:epoch,3,:loss,1.8713183968407767)

### SAMPLE RUN 4ce58d1+ Fri Sep 16 12:24:00 EEST 2016
### Transposed everything so getindex does not need to copy

# charlm.jl (c) Emre Yolcu, Deniz Yuret, 2016. Character level language model based on http://karpathy.github.io/2015/05/21/rnn-effectiveness.
# opts=(:keepstate,false)(:lr,4.0)(:atype,"KnetArray{Float32}")(:winit,0.3)(:savefile,nothing)(:loadfile,nothing)(:generate,0)(:bestfile,nothing)(:gclip,3.0)(:embedding,256)(:hidden,256)(:epochs,3)(:decay,0.9)(:gcheck,0)(:seqlength,100)(:seed,42)(:batchsize,128)(:datafiles,Any["data/10.txt"])(:fast,true)
# INFO: Chars read: [("10.txt",425808)]
# INFO: 87 unique chars.
#   1.388652 seconds (1.73 M allocations: 196.686 MB, 2.04% gc time)
# (:epoch,0,:loss,6.1075509900258)
#   3.940298 seconds (6.06 M allocations: 373.166 MB, 2.06% gc time)
#   3.935995 seconds (6.06 M allocations: 373.244 MB, 2.07% gc time)
#   3.934983 seconds (6.06 M allocations: 373.245 MB, 2.08% gc time)
#   1.390374 seconds (1.73 M allocations: 196.820 MB, 2.11% gc time)
# (:epoch,3,:loss,1.860654126576015)

### SAMPLE RUN 31136d5+ Wed Sep 14 17:51:44 EEST 2016: using vcat(x,h) and vcat(w...)
### optimized learning parameters: winit=0.3, lr=4.0, gclip=3.0

# charlm.jl (c) Emre Yolcu, Deniz Yuret, 2016. Character level language model based on http://karpathy.github.io/2015/05/21/rnn-effectiveness.
# opts=(:keepstate,false)(:lr,4.0)(:atype,"KnetArray{Float32}")(:winit,0.3)(:savefile,nothing)(:loadfile,nothing)(:generate,0)(:bestfile,nothing)(:gclip,3.0)(:embedding,256)(:hidden,256)(:epochs,3)(:decay,0.9)(:gcheck,0)(:seqlength,100)(:seed,42)(:batchsize,128)(:datafiles,Any["10.txt"])(:fast,true)
# INFO: Chars read: [("10.txt",425808)]
# INFO: 87 unique chars.
#   1.596959 seconds (1.79 M allocations: 197.432 MB, 2.56% gc time)
# (:epoch,0,:loss,5.541976199042528)
#   4.421566 seconds (6.23 M allocations: 375.843 MB, 2.54% gc time)
#   4.418540 seconds (6.25 M allocations: 376.058 MB, 2.51% gc time)
#   4.402317 seconds (6.26 M allocations: 376.297 MB, 2.66% gc time)
#   1.594677 seconds (1.81 M allocations: 197.737 MB, 2.64% gc time)
# (:epoch,3,:loss,1.8484550957572192)

### SAMPLE RUN 80503e7+ Wed Sep 14 17:35:36 EEST 2016: using vcat(x,h)
#
# charlm.jl (c) Emre Yolcu, Deniz Yuret, 2016. Character level language model based on http://karpathy.github.io/2015/05/21/rnn-effectiveness.
# opts=(:keepstate,false)(:lr,1.0)(:atype,"KnetArray{Float32}")(:savefile,nothing)(:loadfile,nothing)(:generate,0)(:bestfile,nothing)(:embedding,256)(:gclip,5.0)(:hidden,256)(:epochs,3)(:decay,0.9)(:gcheck,0)(:seqlength,100)(:seed,42)(:batchsize,128)(:datafiles,Any["10.txt"])(:fast,true)
# INFO: Chars read: [("10.txt",425808)]
# INFO: 87 unique chars.
#   1.930180 seconds (1.95 M allocations: 213.741 MB, 1.82% gc time)
# (:epoch,0,:loss,4.462641664662756)
#   4.968101 seconds (7.47 M allocations: 454.259 MB, 2.24% gc time)
#   4.963733 seconds (7.47 M allocations: 454.363 MB, 2.26% gc time)
#   4.967413 seconds (7.45 M allocations: 454.024 MB, 2.14% gc time)
#   1.945658 seconds (1.98 M allocations: 214.183 MB, 2.02% gc time)
# (:epoch,3,:loss,3.2389672966290237)

### SAMPLE RUN 65f57ff+ Wed Sep 14 10:02:30 EEST 2016: separate x, h, w, b
#
# charlm.jl (c) Emre Yolcu, Deniz Yuret, 2016. Character level language model based on http://karpathy.github.io/2015/05/21/rnn-effectiveness.
# opts=(:keepstate,false)(:lr,1.0)(:atype,"KnetArray{Float32}")(:savefile,nothing)(:loadfile,nothing)(:generate,0)(:bestfile,nothing)(:embedding,256)(:gclip,5.0)(:hidden,256)(:epochs,3)(:decay,0.9)(:gcheck,0)(:seqlength,100)(:seed,42)(:batchsize,128)(:datafiles,Any["10.txt"])(:fast,true)
# INFO: Chars read: [("10.txt",425808)]
# INFO: 87 unique chars.
#   2.156358 seconds (2.31 M allocations: 237.913 MB, 2.30% gc time)
# (:epoch,0,:loss,4.465127425659868)
#   6.287736 seconds (9.54 M allocations: 574.703 MB, 2.84% gc time)
#   6.272144 seconds (9.54 M allocations: 574.633 MB, 2.80% gc time)
#   6.277462 seconds (9.54 M allocations: 574.637 MB, 2.86% gc time)
#   2.165516 seconds (2.34 M allocations: 238.323 MB, 2.56% gc time)
# (:epoch,3,:loss,3.226540256084356)

### SAMPLE OUTPUT (with head -10000 100.txt): first version
# julia> CharLM.main("--gpu --data 10.txt")
# opts=(:lr,1.0)(:savefile,nothing)(:loadfile,nothing)(:dropout,0.0)(:generate,0)(:bestfile,nothing)(:embedding,256)(:gclip,5.0)(:hidden,256)(:epochs,10)(:nlayer,1)(:decay,0.9)(:gpu,true)(:seqlength,100)(:seed,42)(:batchsize,128)(:datafiles,Any["10.txt"])
# INFO: Chars read: [("10.txt",425808)]
# INFO: 87 unique chars.
# (0,4.465127425659868)
#   2.182693 seconds (2.36 M allocations: 240.394 MB, 1.74% gc time)
#   7.861079 seconds (10.12 M allocations: 601.311 MB, 1.84% gc time)
# (1,3.3244698245543285)
#   2.159062 seconds (2.35 M allocations: 239.217 MB, 1.80% gc time)
#   6.200085 seconds (9.55 M allocations: 575.573 MB, 2.28% gc time)
# (2,3.24593908969621)
#   2.352389 seconds (2.34 M allocations: 239.381 MB, 1.51% gc time)
#   6.211946 seconds (9.55 M allocations: 575.568 MB, 2.21% gc time)
# (3,3.226540256084356)


text:
using StatsBase
using RDatasets
a = [123,4234,23423,1231231,1432432423,1341413413]
geomean(a)
harmmean(a)
trimmean(a,0.1)

wv = rand(6)
mean(a, weights(wv))

max_stars_repo_path: Chapter 3/4.jl

max_stars_repo_name: PacktPublishing/Julia-for-data-science

max_issues_repo_path: Chapter 3/4.jl

max_issues_repo_name: PacktPublishing/Julia-for-data-science

max_forks_repo_path: Chapter 3/4.jl

max_forks_repo_name: PacktPublishing/Julia-for-data-science


text:
using Documenter
using GeometricProblems

makedocs(
    sitename = "GeometricProblems.jl",
    format = Documenter.HTML(prettyurls = get(ENV, "CI", nothing) == "true"),
    pages = ["Home" => "index.md",
             "Diagnostics"               => "diagnostics.md",
             "Exponential Growth"        => "exponential_growth.md",
             "Hénon-Heiles System"       => "henon_heiles.md",
             "Kepler Problem"            => "kepler_problem.md",
             "Lorenz Attractor"          => "lorenz_attractor.md",
             "Lotka-Volterra 2d"         => "lotka_volterra_2d.md",
             "Lotka-Volterra 3d"         => "lotka_volterra_3d.md",
             "Lotka-Volterra 4d"         => "lotka_volterra_4d.md",
             "Massless Charged Particle" => "massless_charged_particle.md",
             "Harmonic Oscillator"       => "harmonic_oscillator.md",
             "Nonlinear Oscillators"     => "nonlinear_oscillators.md",
             "Mathematical Pendulum"     => "pendulum.md",
             "Double Pendulum"           => "double_pendulum.md",
             "Point Vortices"            => "point_vortices.md",
             "Inner Solar System"        => "inner_solar_system.md",
             "Outer Solar System"        => "outer_solar_system.md",
            ]
)

deploydocs(
    repo   = "github.com/JuliaGNI/GeometricProblems.jl",
)