refactor icethermo0 and icethermoh

[fjansson]

I have refactored the routines icethermo0 and icethermoh in modthermodynamics.f90:

• a separate function for determining qsatur - this was inlined in the Newton's method code and hard to read
• added some comments
• added comments and a reference in modglobal.f90, where the saturation tables are calculated

Would you like to merge this? The commit is here: CloudResolvingClimateModeling@53c38ca

Testing:
I compared the prognostic fields after 1000s of simulation (cases/example) with and without this change, the sum(difference**2) was < 1e-8 for all of them. Is there some official testing procedure?

This should be possible to merge into the Dales master branch. I didn't make a real pull request, since that would include all previous commits in our branch as well. I can create a new branch and make a real pull request - which dales branch should I base it on then?