How to color system of PDEs with two population

Dear Matlab Team,

I have a system of PDEs with two equations (two population, V, W) written in polar coordinate. My diffusion is nonhomogeneous and I get the centripetal pattern from each population. For example from 0 degree to 30 degree is V and from 30 to 60 is W and this process repeates. Now, I need to color V as a red and W as a blue in a circle and I use this code, but since the green channel is zero, it is black between grid points.

Could you help me how to finx this to have a smoth diffusion with no black distance?

 [Theta, R] = meshgrid(theta, r);  % [51, 60]
       [X, Y] = pol2cart(Theta, R);      % [51, 60]
       Vnorm = V / (max(V(:)) + 1e-8);
       Wnorm = W / (max(W(:)) + 1e-8);
      RGB = zeros(Nr+1, Ntheta, 3);     % [r, theta, color channel]
      RGB(:,:,1) = Wnorm;   % 
      RGB(:,:,2) = 0;       % 
      RGB(:,:,3) = Vnorm;   % 

     clf;
     surf(X, Y, zeros(size(X)), RGB, 'EdgeColor', 'none');
     view(2);
     axis equal off;
     title(['V (Blue) , W (Red), t = ', num2str(t*dt)]);
     colorbar;
     drawnow;

Thank you !!

[William Rose]

@Neda,

It would help to have runnable code, so we can actually see the problem. Your posted code is un-runnable since various scalars, vectors, and arrays are not provided. Please provide the simplest possible runnable code and run it in the Matlab Answers window. For example, you could save the necessary variables in a .mat file, attach it, then load it at the start of the script posted on Matlab Answers.

[William Rose]

@Neda,

You said "For example from 0 degree to 30 degree is V and from 30 to 60 is W and this process repeates." I'm not sure I understand what you mean. Here is an example of V and W, where V is zero from 30 to 60 deg, 90 to 120, etc.; and W is zero from 0 to 30, 60 to 90, etc. Is that what you meant?

r=(0:0.2:10)'; Nr=length(r);
theta=0:pi/36:2*pi; Nt=length(theta);
Vr=tukeywin(Nr); Vt=max(sin(6*theta),0);
Wr=tukeywin(Nr); Wt=max(-sin(6*theta),0);
V=Vr*Vt; 
W=Wr*Wt;
% Plot Vr, Wr, Vt, Wt
figure
subplot(211), plot(r,Vr,'-rx',r,Wr,'-bo')
xlabel('Radius'); legend('V(r)','W(r)')
subplot(212), plot(theta,Vt,'-rx',theta,Wt,'-bo')
xlabel('Theta'); legend('V(theta)','W(theta)')

By the way, it looks like your posted code assigns V to blue and W to red, which is the opposite of what you said you wanted. No big deal.

Answer 1

[William Rose]

@Neda,

Using my code from the comment above, we can run the rest of your code:

r=(0:0.2:10)'; Nr=length(r);
theta=0:pi/36:2*pi; Ntheta=length(theta);
Vr=tukeywin(Nr); Vt=max(sin(6*theta),0);
Wr=tukeywin(Nr); Wt=max(-sin(6*theta),0);
V=Vr*Vt; 
W=Wr*Wt;
[Theta, R] = meshgrid(theta, r);  % [51, 60]
[X, Y] = pol2cart(Theta, R);      % [51, 60]
Vnorm = V / (max(V(:)) + 1e-8);
Wnorm = W / (max(W(:)) + 1e-8);
RGB = zeros(Nr, Ntheta, 3);     % [r, theta, color channel]
RGB(:,:,1) = Wnorm;   % 
RGB(:,:,2) = 0;       % 
RGB(:,:,3) = Vnorm;   % 
clf;
surf(X, Y, zeros(size(X)), RGB, 'EdgeColor', 'none');
view(2);
axis equal off;
title('V (Blue) , W (Red)');

●

I think this plot is doing what you want. It is consistent with the plots of V verus r, V vs. theta, W vs. r, W vs. theta, in my earlier post. There is not black between the grid points. There is black around theta=0, 30, 60... because both V and W are zero at those angles. I deleted the colorbar, since it was misleading. You would need two colorbars -one for V and one for W - and I don't know how to do that, or if it can be done.

[Neda]

Thank you so much @William Rose

Is it possible to write this code without this part:

RGB = zeros(Nr, Ntheta, 3); % [r, theta, color channel]

RGB(:,:,1) = Wnorm; %

RGB(:,:,2) = 0; %

RGB(:,:,3) = Vnorm; %

[Walter Roberson]

RGB = cat(3, Wnorm, zeros(size(Wnorm)), Vnorm);

[William Rose]

@Neda, you're welcome.

@Walter Roberson, Thanks. Nice.

[William Rose]

@Neda,

Now I'll have a little fun. In the script below I do the original idea for W (called W1 below), where it repeats every 60 degrees, and I try it wher W repeats every 45 degrees (W2). This causes magenta spokes to appear in the plot, where W2 and V are both maximum at about the same angles. V is unchanged.

r=(0:0.2:10)'; Nr=length(r);
theta=0:pi/72:2*pi; Ntheta=length(theta);
Vr=tukeywin(Nr); Vt=max(sin(6*theta),0);
Wr=tukeywin(Nr); Wt1=max(-sin(6*theta),0); Wt2=max(-sin(8*theta),0);
V=Vr*Vt; 
W1=Wr*Wt1;
W2=Wr*Wt2;
[Theta, R] = meshgrid(theta, r);  % [51, 60]
[X, Y] = pol2cart(Theta, R);      % [51, 60]
Vnorm = V / (max(V(:)) + 1e-8);
Wnorm1 = W1 / (max(W1(:)) + 1e-8);
Wnorm2 = W2 / (max(W2(:)) + 1e-8);
RGB1 = cat(3, Wnorm1, zeros(size(Vnorm)), Vnorm);
RGB2 = cat(3, Wnorm2, zeros(size(Vnorm)), Vnorm);
% Plot Vr, Wr, Vt, Wt
figure
subplot(211), plot(r,Vr,'-rx',r,Wr,'-bo')
xlabel('Radius'); legend('V(r)','W(r)')
subplot(212), plot(theta,Vt,'-rx',theta,Wt1,'-bo',theta,Wt2,'-b.')
xlabel('Theta'); legend('V(theta)','W1(theta)','W2(theta)')

% Plot V, W1 (polar)
figure
subplot(121), surf(X, Y, zeros(size(X)), RGB1, 'EdgeColor', 'none');
view(2);
axis equal off;
title('V (Blue) , W1 (Red)');
% Plot V, W2 (polar)
subplot(122), surf(X, Y, zeros(size(X)), RGB2, 'EdgeColor', 'none');
view(2);
axis equal off;
title('V (Blue) , W2 (Red)');