Plot velocity profile of microfluid flow

Hi,

My project is in microfluidics and I am trying to rewrite the following velocity profile into matlab so that I can control XY stepper motors accordingly. First of all, I want to turn the below equations into code that is then plotted on a graph of velocity against position.

This velocity profile depicts the velocity of a microfluid moving through a converging/diverging microchannel, where xdim represents dimensionless position(x) in the channel, u represents velocity and f1, f1, n1 and n2 are dimensionless parameters.

Below is the code I have written so far, although I am new to matlab so I am unsure where I have gone wrong here and I don't know if using a for loop if the correct thing to do.

for idx = 1:6
    if xdim <= -n1 - n2/2
       u = Uu
       plot(xdim,u,"-") 
    elseif xdim <= -n1 - n2/2
        u = f2*(xdim + n1 + n2/2)^2 + Uu
        plot(xdim,u,"-") 
    elseif -n1 + n2/2 <= xdim <= 0
        u = f1*xdim + Uc
        plot(xdim,u,"-") 
    elseif 0 <= xdim <= n1 - n2/2
        u = f1*xdim + Uc
        plot(xdim,u,"-") 
    elseif 0 <= xdim <= n1-m2/2
        u = f2*(xdim-n1-n2/2)^2 + Uu
        plot(xdim,u,"-") 
    elseif xdim >= n1 + n2/2
        u = Uu
        plot(xdim,u,"-")
    end
        hold on
end

Any help is greatly appreciated!

Thank you.

[Lewis Kirkwood]

Aplogies, the typed equations were incorrect. Please find the correct equations below. However, I do not receive the correct plot when running the code

for idx = 1:6
    if xdim <= -n1 - n2/2
       u = Uu
       plot(xdim,u,"-") 
    elseif -n1 - n2/2 < xdim < -n1 - n2/2
        u = f2*(xdim + n1 + n2/2)^2 + Uu
        plot(xdim,u,"-") 
    elseif -n1 + n2/2 <= xdim <= 0
        u = f1*xdim + Uc
        plot(xdim,u,"-") 
    elseif 0 <= xdim <= n1 - n2/2
        u = -f1*xdim + Uc
        plot(xdim,u,"-") 
    elseif n1 - n2/2 < xdim < n1-n2/2
        u = f2*(xdim-n1-n2/2)^2 + Uu
        plot(xdim,u,"-") 
    elseif xdim >= n1 + n2/2
        u = Uu
        plot(xdim,u,"-")
    end
        hold on

end

Answer 1

[Jon]

Your use of the if elseif to determine what the appropriate function to evaluate is good, but you have a conceptual error on what you are looping through.

Basically you need to first define a vector xdim of the values of the dimensionless distance you are interested in.Then you should loop through each element of this vector, determine the range it is in, and the corresponding function for calculating the velocity, u, and then calculate it.

Finally once you have all of the u values calculated you plot u vs xdim.

So something like this:

% you'll need to assign right hand side for all of these variables
f1 = 
f2 =
Uu =
Uc =
n1 = 
n2 =
xdimMin = % smallest value of xdim
xdimMax = % largest value of xdim
numPoints = % number of xdim values that you want to evaluate between xdimMin and xdimMax

% define your xdim values
xdim = linspace(xdimMin,xdimMax,numPoints) % numPoints equally space values

% preallocate a vector to hold your velocities (not necessary but good practice)
u = zeros(numPoints,1)

% loop through xdim values calculating velocities using appropriate function for the range of xdim
for idx = 1:numPoints
    if xdim(idx) <= -n1 - n2/2
       u(idx) = Uu

    elseif xdim(idx)< -n1 + n2/2
        u(idx) = f2*(xdim + n1 + n2/2)^2 + Uu

    elseif xdim(idx)<= 0
        u(idx) = f1*xdim + Uc

    elseif xdim(idx) <= n1 - n2/2
        u(idx) = -f1*xdim + Uc

    elseif xdim(idx) < n1 + n2/2
        u(idx) = f2*(xdim-n1-n2/2)^2 + Uu

    elseif xdim(idx)>= n1 + n2/2
        u(idx) = Uu

    end      
end

% plot results
plot(xdim,u)
xlabel('dimensionless position xdim')
ylabel('dimensionless velocity u')

Alternatively you could just have a line of code for each range of xdim and assign all of the u values for the x's in each range and then plot so something like

% assume you've already assigned n1,n2,f1,...
% define your xdim values
xdim = linspace(xdimMin,xdimMax,numPoints); % numPoints equally space values

% preallocate a vector to hold your velocities (not necessary but good practice)
u = zeros(numPoints,1);

% calculate velocity using appropriate function for range of xdim
% note use of logical indexing on left hand side
u(xdim <= -n1 - n2/2) = Uu;
u(-n1 - n2/2 < xdim & xdim< -n1 - n2/2) = f2*(xdim + n1 + n2/2)^2 + Uu;
u( -n1 + n2/2 <= xdim & xdim <= 0) = f1*xdim + Uc;
u( 0 <= xdim & xdim<= n1 - n2/2) = -f1*xdim + Uc;
u(n1 - n2/2 < xdim & xdim < n1-n2/2) = f2*(xdim-n1-n2/2)^2 + Uu;
u(xdim >= n1 + n2/2) = Uu;

% plot results
plot(xdim,u)
xlabel('dimensionless position xdim')
ylabel('dimensionless velocity u')

p.s. Better check the ranges I have carefully to make sure there aren't some typos in there

[Lewis Kirkwood]

Thank you Jon, I can now plot the movement I desire in Matlab.

Following on, I want to command my stage controller (Prior ProScan III) to move the stage according to this velocity profile in the X-axis only. I can use the 'Instrument Control' app to connect to my controller via a serial USB port but from there I am stuck as to how to get it to follow the code.

Do you have any advice on this?

[Jon]

Glad to hear that this helped you.

I'm sorry I don't have experience with the "Instrument Contro" app. Since it sounds like you are now moving on to a new problem, I suggest trying to clearly frame this follow up question and and post it as a new self contained question in MATLAB answers