When modelling the melting of snow on sea-ice, how can I create a vector, for example the density of both snow and ice, respectively, so that I can run it in a loop?

I am currently writing a code to model the evolution of melt-ponds on sea-ice. For simplicity, I have divided the snow and ice into 5 layers each. I am using a diffusion equation. My aim is to model the melt within a loop, so I believe I need to create vectors for constants like specific heat capacity and density for the two given mediums. The density for snow = 200 and ice = 919 (kg m-3). So I can then calculate this in the model

(dt/(rho.*c.*dz))*Q

Thank you in advance.

Answer 1

[dpb]

...show me how a vector might be constructed for that? (density: snow=200, ice=919, water=1000)

Basically, it's an exercise in indexing by array in Matlab...

>> rhos=[200;1000;919];  % assign densities to states of 1:3--> snow, water ice
>> state_v=repmat([1:3].',1,5).'; state_v=state_v(:);  % just make up something
>> rho=rhos(v)
rho =
       200
       200
       200
       200
       200
      1000
      1000
      1000
      1000
      1000
       919
       919
       919
       919
       919
>>

So you have the one array of the basic property by material and then use that in a lookup table to create the actual property for each layer for the computation. Again, you just check for a state transisition for each layer at each time step before compute the new overall vector.

I also showed another manipulation "trick" w/ Matlab of rearranging arrays and vectors to get some other desired shape and/or order. The state I start with here is likely quite unrealistic with all three states of equal numbers but it's just again, an example for illustration only.

You can, of course, do the same thing with only two states and in fact, may be simpler than the logical I started with. That came to me originally as it was the binary case but even so a 1:2 choice instead of T|F would let you use the same idea as here.

This might be a place where a structure could be of use in factoring the model into a given set of parameters for a state but that's getting a little farther afield from your actual query.

Answer 2

[dpb]

Well, it's easy enough to create the vectors as

rho=200*ones(5,1);

say, for snow but I'm presuming you're keeping a moving boundary layer of which layer is ice/snow in the five? Otherwise, you can always multiply any array by a constant.

So, looks like what you'll need is some logical array of

isIce=true(5,1);  % initially all ice

and then keep the density vectors as something like

rho(isIce)=919; rho(~isIce)=200;

For convenience or legibility or just because you could do have

isSnow=~isIce;

with which you can avoid the not expression in the subscript expression for some legibility.

ADDENDUM

And, of course, you can get slightly more sophisticated code-wise at some loss of transparency in what is actually going on and write things like --

rho=200+719*isIce;

and each time step update the isIce logical vector then compute the new density vector.

>> isIce=rand(5,1)>0.5
isIce =
   0
   1
   1
   0
   1
>> rho=200+719*isIce
rho =
 200
 919
 919
 200
 919
>>

[Thomas]

Thank you for your reply. That does make sense. So I have got two separate sets of vectors; five 200s and five 919s for rho(~isIce) and rho(isIce), respectively.

How would I 'combine' these two vectors into one rho vector, with snow on top of ice, i.e. five 200s above five 900s? If that makes sense. Thanks again.

[dpb]

Make the logical vector as long as it needs be and set the values of it based on the criterion of which level is/isn't ice. Same thing as above where I randomized a length of five but the initial state would be

isIce=[false(5,1);true(5,1)];
rho=200+719*isIce;

You can, of course, also use a named variable for the two densities instead of "magic numbers". That could be something like --

rhoIce=919;
rhoSno=200;
rho=rhoSno + (rhoIce-rhoSno)*isIce;

[dpb]

Actually seems like you'd need three states instead of only two as surely there's molten water in the mix at times, too, isn't there? If so, the same idea works except instead of just a binary logical you'll build a state vector of three values and fix up a similar interpolating/computing/lookup scheme to set the values as I used above for the binary case.

[Thomas]

To your first reply: I have converted what I need into logic vectors - it seems to work well.

Your second reply: Yes, ideally, that would be the plan. So in my (very) simplified model, when a layer of snow melts I'm going to treat it as the water percolates into the layer below and runs off. So I don't want there to be a layer of water sitting on top of the snow. When it comes to the ice, then yes, I would like to represent water overlying.

Could you possibly show me how a vector might be constructed for that? (density: snow=200, ice=919, water=1000)

Thanks, that would be really appreciated.