#!/usr/bin/ruby
# File: lab5_v3.rb
# => Resource: Input5.txt

#Initialize Arrays:
M = Array.new(11) {0.0}
EM = Array.new(11) {0.0}
D = Array.new(11) {0.0}
W = Array.new(11) {0.0}
B = Array.new(11) {0.0}
ED = Array.new(11) {0.0}
MO = Array.new(11) {0.0} 
C = Array.new(11) {0.0}
EC = Array.new(11) {0.0}

F = Array.new(11) {0.0}
EF = Array.new(11) {0.0}
P = Array.new(11) {0.0}
G = Array.new(11) {0.0}
EG = Array.new(11) {0.0}
I = Array.new(11) {0.0}
EI = Array.new(11) {0.0}
EMO = Array.new(11) {0.0}

J = Array.new(11) {0.0}
EJ = Array.new(11) {0.0}
EP = Array.new(11) {0.0} 
Z = Array.new(11) {0.0}
MS = Array.new(11) {0.0}
EMS = Array.new(11) {0.0}
X = Array.new(11) {0.0}
V1 = Array.new(11) {0.0}

L = Array.new(11) {0.0}
EL = Array.new(11) {0.0}
V10 = Array.new(11) {0.0}
EV10 = Array.new(11) {0.0}
SV1 = Array.new(11) {0.0}
V2 = Array.new(11) {0.0}
MP = Array.new(11) {0.0}
EMP = Array.new(11) {0.0}

SV2 = Array.new(11) {0.0}
MI = Array.new(11) {0.0}

line = File.readlines("Input5.txt")
name = line[0]
experiment = line[1]
partner = line[2]
date = line[3]
empty_pycnometer = line[4].split(',')
we = empty_pycnometer[0].to_f
ewe = empty_pycnometer[1].to_f

pycometer_water = line[5].split(',')
wo = pycometer_water[0].to_f
ewo = pycometer_water[1].to_f

t = line[6]
density_water = line[7].split(',')
d0 = density_water[0].to_f
ed0 = density_water[1].to_f

nacl = line[8].split(',')
M[1] = nacl[0].to_f
EM[1] = nacl[1].to_f

rum = EM[1]/M[1]
n = line[9].to_i

out = File.open("Output5.txt", "w")
out.printf("NAME: %s",name) 
out.printf("EXPERIMENT: %s",experiment)
out.puts "PARTNER: #{partner}"
out.puts "DATE: #{date}"
out.printf("THE WEIGHT OF THE EMPTY PYCOMETER = %d7.4 +/- %d5.4 grams\n", we, ewe)
#out.puts "  = #{we} +/- #{ewe} grams"
out.puts "THE WEIGHT OF THE PYCNOMETER FILLED WITH H2O"
out.puts "  = #{wo} +/- #{ewo} grams"
out.puts "THE TEMPERATURE OF THE BATH"
out.puts "  = #{t} degC"
out.puts "THE DENSITY OF WATER AT #{t} degC"
out.puts "  = #{d0} +/-  #{ed0} grams"
out.puts "THE INITIAL CONCENTRATION OF NACL (MOLARITY)"
out.puts "  = #{M[1]} +/- #{EM[1]} mol/dm3"
#for i in 1..n
out.puts "#{n} CONCENTRATIONS WERE RUN"

vp = (wo-we)/d0  #vp = VOLUME OF THE PYCNOMETER
a = wo-we
ea = (ewo**2 + ewe**2)**0.5
evp = vp * ((ea/a) **2 + (ed0/d0) **2) **0.5
#        ERROR IN THE VOLUME OF THE PYCNOMETER

ewp = line[10].to_f
out.print "THE UNCERTAINTY IN THE WEIGHT OF THE PYCNOMETER FILLED "
out.puts "WITH THE VARIOUS NACL SOLUTIONS = +/- #{ewp} grams\n\n"

for i in 2..n
  M[i] = M[i-1]/2.0
  EM[i] = rum * M[i]
#puts "M[#{i}] = #{M[i]} and EM[#{i}] = #{EM[i]}"
end

current = 11 # line containing mass of first concentration; rest follow
for i in 1..n
  out.puts
  out.puts "RUN NO. #{i}"
  out.print "THE CONCENTRATION (MOLARITY) FOR \##{i}"
  out.puts " = #{M[i]} +/- #{EM[i]} mol/dm3"
  W[i] = line[current].to_f
  out.puts "THE WEIGHT OF THE FILLED PYCNOMETER = #{W[i]} grams"
  current = current +1
end

for i in 1..n
  out.puts
  #D[i] = W[i]-we/vp
  # =DENSITIES
  
  B[i] = W[i] - we
  eb = (ewp**2 + ewe**2)**0.5
  D[i] = B[i] / vp
  ED[i] = D[i] * ((eb/B[i])**2 + (evp/vp)**2)**0.5 
  # => ERROR IN THE DENSITY OF THE NACL SOLUTIONS
  # MO[i]=1/((D[i]/M[i]) - (58.45/1000))
  # => MOLALITIES
  C[i]=D[i] / M[i]
  EC[i] = C[i] * ((ED[i]/D[i])**2 + (EM[i] / M[i])**2)**0.5
  F[i] = C[i] - 58.45/1000
  EF[i] = EC[i]
  MO[i] = 1.0 / F[i]
  EMO[i] = MO[i] * EF[i] / F[i]
  # => ERROR IN THE MOLALITIES OF THE NACL SOLUTIONS
  #
  #P(I)=(1/D(I))*(58.45-(1000/MO(I))*(W(I)-WO)/(WO-WE))
  # => PHI
  #
  G[i] = W[i] - wo
  EG[i] = (ewp**2 + ewo**2)**0.5
  h = wo - we
  eh = (ewo**2 + ewe**2)**0.5
  I[i] = 1000 * G[i] / (MO[i] * h)
  EI[i] = I[i] * ((EG[i] / G[i])**2 + (EMO[i] / MO[i])**2 + (eh / h)**2)**0.5
  J[i] = 58.45 - I[i]
  EJ[i] = EI[i]
  P[i] = J[i] / D[i]
  EP[i] = P[i] * ((EJ[i] / J[i])**2 + (ED[i] / D[i])**2)**0.5
  # => ERROR IN PHI
  #
  # I[i]=(1000/MO[i])*(W[i]-wo)/(wo-we)        (1000/m)(w-wo)/(wo-we)
  # => FUNCTION (TO BE PLOTTED)
  #
  # EI[i]= ERROR IN THIS FUNCTION TO BE PLOTTED
  #
end

for p in 1..n
  Z[p] = MO[p]**0.5   #Z[p] = SQUARE ROOT OF MOLALITY
  MS[p] = Z[p]
  EMS[p] = 0.5 * Z[p] * EMO[p] / MO[p]
  # => ERROR IN THE SQUARE ROOT OF MOLALITY
end

sumx = 0
sumy = 0
sumy2 = 0
sumx2 = 0
sumxy = 0

for ii in 1..n
  sumx = sumx + Z[ii]
  sumy = sumy + P[ii]
  sumy2 = sumy2 + P[ii]**2
  sumx2 = sumx2 + Z[ii]**2
  sumxy = sumxy + P[ii]*Z[ii]
end

xmean = sumx / n
ymean = sumy / n
yint = (sumx2*sumy - sumx*sumxy) / (n*sumx2 - sumx**2)      # y-int = PHI
slope = (sumy - n*yint) / sumx                              #slope
qq = n * sumxy - sumx*sumy
pp = ((n*sumx2 - (sumx**2)) * (n*sumy2 - (sumy**2)))**0.5
r = qq/pp           # (CORRELATION COEFFICIENT)
r2 = (r**2) * 100   # (PERCENT OF FIT)

sni = 0 
sigma = 0
for k in 1..n
  ycalc = slope * Z[k] + yint
  sni = (ycalc - P[k])
  sni = sni.abs
  sigma = sigma + sni**2
end

sy = (sigma / (n - 2))**0.5                       # sy = ERROR IN Y
sm = ((n / (n * sumx2 - sumx**2))**0.5) * sy      # sm = ERROR IN SLOPE
sb = ((sumx2 / (n * sumx2 - sumx**2))**0.5) * sy  # sb = ERROR IN Y-INT
X[1] = 0                                          # X[I] = MOLALITY

for i in 1..6
  # V1(I)=(18.016/D0)-((X(I)/55.51)*(X(I)^0.5)*(SLOPE/2))
  # => PARTIAL MOLAR VOLUME OF V1
  #
  L[i] = X[i]**1.5 * slope / (55.506 * 2)
  EL[i] = X[i]**1.5 * sm / (55.506 * 2)
  v1O = 18.016 / 0.997
  ev1O = v1O * ((0.001 / 18.016)**2 + (0.0001 / 0.997)**2)**0.5
  V1[i] = v1O - L[i]
  SV1[i] = (ev1O**2 + EL[i]**2)**0.5
  # => ERROR IN PARTIAL MOLAR VOLUME OF V1
  #
  # V2(I) = YINT + 1.5*SLOPE*X(I)^0.5
  # => PARTIAL MOLAR VOLUME OF V2
  #
  MP[i] = 1.5 * slope * X[i]**0.5
  EMP[i] = 1.5 * sm * X[i]**0.5
  V2[i] = yint + MP[i]
  SV2[i] = (sb**2 + EMP[i]**2)**0.5
  # => ERROR IN PARTIAL MOLAR VOLUME OF V2
  #
  X[i + 1] = X[i] + 0.5
end

5.times{out.puts}

out.puts "PLOT OF PHI VS. SQUARE ROOT OF MOLALITY"
out.puts "===========     ======================="
out.puts
out.puts
out.puts "PHI\t\t\t\t\t\tSQUARE ROOT OF MOLALITY"
out.puts
for i in 1..n
  out.puts"#{P[i]}  +/-  #{EP[i]}\t #{MS[i]}  +/-  #{EMS[i]}"
end

out.puts
out.puts "\tSLOPE"
out.puts "\t====="
out.puts "\t  = #{slope} +/- #{sm}"

out.puts
out.puts "\tPHI (O)"
out.puts "\t======="
out.puts "\t  = #{yint} +/- #{sb}"

out.puts
out.puts "\tCORRELATION COEFFICIENT"
out.puts "\t======================="
out.puts "\t  = #{r}"

out.puts
out.puts "\tPERCENT OF FIT"
out.puts "\t=============="
out.puts "\t  = #{r2}%"

3.times{out.puts}

out.puts "PLOT OF V2 VS. MOLALITY"
out.puts "==========     ========"
out.puts
out.puts "V2 \t\t\t\t\t\tMOLALITY"
out.puts "=====================\t\t\t\t========================="
for i in 1..6
  out.print "#{V2[i]} +/- #{SV2[i]}"
  out.printf("%13s", "#{X[i]}\n")
end

3.times{out.puts}

out.puts "PLOT OF V1 VS. MOLALITY"
out.puts "==========     ========"
out.puts
out.puts "V1 \t\t\t\t\t\tMOLALITY"
out.puts "=====================\t\t\t\t========================="
for i in 1..6
  out.print "#{V1[i]} +/- #{SV1[i]}"
  out.printf("%13s", "#{X[i]}\n")
end

3.times{out.puts}

out.puts "\tRUN NO.\t DENSITY"
out.puts "\t=======\t ================================="
for i in 1..n
  out.puts "\t #{i}\t #{D[i]}    +/- #{ED[i]}"
end
out.puts
out.puts

out.puts "\tRUN NO.\t MOLARITY"
out.puts "\t=======\t ================================="
for i in 1..n
  out.puts "\t #{i}\t #{M[i]}    +/- #{EM[i]}"
end
out.puts
out.puts

out.puts "\tRUN NO.\t MOLALITY"
out.puts "\t=======\t ================================="
for i in 1..n
  out.puts "\t #{i}\t #{MO[i]}    +/- #{EMO[i]}"
end
out.puts
out.puts

out.puts "\tRUN NO.\t PHI"
out.puts "\t=======\t ================================="
for i in 1..n
  out.puts "\t #{i}\t #{P[i]}    +/- #{EP[i]}"
end
out.puts
out.puts

out.puts "\tRUN NO.\t (1000/M) (W-WO) / (WO-WE)"
out.puts "\t=======\t ================================="
for i in 1..n
  out.puts "\t #{i}\t #{I[i]}    +/- #{EI[i]}"
end
out.puts
out.puts

out.puts "VOLUME OF PYCNOMETER = #{vp} +/- #{evp} CM3"
out.puts "====== == =========="
out.puts
out.puts
out.puts "PHI (O) = #{yint} +/- #{sb} CM3"
out.puts "======="
out.puts
out.puts
out.puts "DPHI/DM(1/2) = #{slope} +/- #{sm} CM3/MOLALITY(1/2)"
out.puts "============"
 
exit