3 "for" loops optimization or changing for global maximum finding

Hi,

I have a function (function1), which has three variables (i,j,k) which I need to change from -20 to 20 in step 1, and later with function2, I have to find maximum DB number and place (ci, cj, ck), where DB is the highest (it can be everyhere, depending on initial data file).

This code works good, but the main disadvantage of it is time. Calculations take about 3 hours (a lot of combinations ofc.). I need somehow to optimise my code. Can someone explain how can I find global maximum of DB faster, than doing 3 "for" loops?

I have heard about some algorithms which works, but I don't know how and how to use it in matlab.

Thanks.

max = 0;
 ci = 0;
 cj = 0;
 ck = 0;
%  kr1 = 0;
%  kg2 = 0;
%  kb3 = 0;

for i = -20:1:20
    for j = -20:1:20
        for k = -20:1:20
            [spd, kr, kg, kb] = function1( i, j, k, A, B, C );   
            DB = function2( mspd4, mspd5, Xn, Yn, Zn, U, range, D, E, F, spd );
            if(max<=DB) 
                max = DB;
                ci = i; cj = j; ck = k;
               % kr1 = kr; kg2 = kg; kb3 = kb; 
            end
        end
    end
end

[Michael Croucher]

Can you share function1 and function2 to make this an example that's runnable please?

[Domantas B]

Function1:

function [spd,kr,kg,kb] = spektrustumdymas( kieknm,kieknm2,kieknm3, A,B,C)

range1 = 451;
range2 = 680;
rangered = range1:range2;

unitconversion = (1e-6/.0001);


 spektras1 = mspd1; %initial file 1
 spektras2 = mspd2; %iniitial file 2
 spektras3 = mspd3; %initial file 3
 
% for i = 1:471 
   %  spektras1iki471(i,:) = spektras1(i,:); A
   %  spektras2iki471(i,:) = spektras2(i,:); B
 %    spektras3iki471(i,:) = spektras3(i,:); C
%end

xconst = 0.3608;
yconst = 0.3635;

rangespektras = A(92:1:321,:);
outspektras1 = A(kieknm+92:1:kieknm+321,:); 
outspektras2 = B(kieknm2+92:1:kieknm2+321,:); 
outspektras3 = C(kieknm3+92:1:kieknm3+321,:); 
outspektras1nm = [outspektras1(:,1)-kieknm, outspektras1(:,2)];
outspektras2nm = [outspektras2(:,1)-kieknm2, outspektras2(:,2)];
outspektras3nm = [outspektras3(:,1)-kieknm3, outspektras3(:,2)];

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% PART FOR XYZ CALCULATIONS:
XYZ2degred = gettristimulus2deg(interpolatespd([outspektras1nm(:,1) outspektras1nm(:,2).*unitconversion],rangered),rangered);
xyz2degred = getxyz(XYZ2degred);

raudonaX = XYZ2degred(1,1);
raudonaY = XYZ2degred(1,2);
raudonaZ = XYZ2degred(1,3);

raudXYZsuma = raudonaX+raudonaY+raudonaZ;

XYZ2deggreen = gettristimulus2deg(interpolatespd([outspektras2nm(:,1) outspektras2nm(:,2).*unitconversion],rangered),rangered);
xyz2deggreen = getxyz(XYZ2deggreen); 

zaliaX = XYZ2deggreen(1,1);
zaliaY = XYZ2deggreen(1,2);
zaliaZ = XYZ2deggreen(1,3);

zaliaXYZsuma = zaliaX+zaliaY+zaliaZ;

XYZ2degblue = gettristimulus2deg(interpolatespd([outspektras3nm(:,1) outspektras3nm(:,2).*unitconversion],rangered),rangered);
xyz2degblue = getxyz(XYZ2degblue);

melynaX = XYZ2degblue(1,1);
melynaY = XYZ2degblue(1,2);
melynaZ = XYZ2degblue(1,3);

melynaXYZsuma = melynaX+melynaY+melynaZ;

raudonaXnormuotas = raudonaX/(raudXYZsuma); 
raudonaYnormuotas = raudonaY/(raudXYZsuma); 
raudonaZnormuotas = raudonaZ/(raudXYZsuma);

zaliaXnormuotas  = zaliaX/(zaliaXYZsuma);
zaliaYnormuotas  = zaliaY/(zaliaXYZsuma);
zaliaZnormuotas  = zaliaZ/(zaliaXYZsuma);

melynaXnormuotas = melynaX/(melynaXYZsuma);
melynaYnormuotas = melynaY/(melynaXYZsuma);
melynaZnormuotas = melynaZ/(melynaXYZsuma);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%% COEFFICIENT Kn calculations:
syms kr kg kb
eqn1 = raudonaXnormuotas*kr + zaliaXnormuotas*kg + melynaXnormuotas*kb == xconst;
eqn2 = raudonaYnormuotas*kr + zaliaYnormuotas*kg + melynaYnormuotas*kb == yconst;
eqn3 = kr + kg + kb == 1;
eqns = [eqn1 eqn2 eqn3];

sol = solve(eqns, [kr, kg, kb]);
kr = double(sol.kr);
kg = double(sol.kg);
kb = double(sol.kb);

spd = [rangespektras(:,1), outspektras1nm(:,2)/(raudonaX+raudonaY+raudonaZ)*kr + outspektras2nm(:,2)/(zaliaX+zaliaY+zaliaZ)*kg + outspektras3nm(:,2)/(melynaX+melynaY+melynaZ)*kb];

[Domantas B]

Function2:

function [ deltabendras ] = dlabbandiniu( mspd4,mspd5,Xn,Yn,Zn,U,range,D,E,F,spd )

sandauga = D(:,2).*spd(:,2); 
naujasspektras = [spd(:,1), sandauga];

sandaugatrwat = E(:,2).*spd(:,2);
naujasspektrastrwat = [spd(:,1), sandaugatrwat];

sandaugacbat = F(:,2).*spd(:,2); 
naujasspektrascbat = [spd(:,1), sandaugacbat];

modelinis = mspd5(92:1:321,:);

modelinistrwat = E(:,2).*modelinis(:,2);
newmodelinistrwat = [modelinis(:,1), modelinistrwat];

modeliniscbat = F(:,2).*modelinis(:,2);
newmodeliniscbat = [modelinis(:,1), modeliniscbat];
%function changing from xyz to lab
 a = labfunkcija(newmodelinistrwat,Xn,Yn,Zn,unitconversion,range); 
 a2 = labfunkcija(newmodeliniscbat,Xn,Yn,Zn,unitconversion,range); 
 a3 = labfunkcija(naujasspektrastrwat,Xn,Yn,Zn,unitconversion,range); 
 a4 = labfunkcija(naujasspektrascbat,Xn,Yn,Zn,unitconversion,range); 
 %%%%
 
 damodelinis = (a2(1,2) - a(1,2));
 dbmodelinis = (a2(1,3) - a(1,3));

 damano = (a4(1,2) - a3(1,2));
 dbmano = (a4(1,3) - a3(1,3));

 deltamodelinis = sqrt((damodelinis)^2 + (dbmodelinis)^2);
 deltamano = sqrt((damano)^2 + (dbmano)^2);

 deltabendras = deltamano - deltamodelinis;

end

[Walter Roberson]

            [spd, kr, kg, kb] = function1( i, j, k, A, B, C );   
            DB = function2( mspd4, mspd5, Xn, Yn, Zn, U, range, D, E, F, spd );

Most of those variables are not defined. A, B, C, mspd4, mspd5, Xn, Yn, Zn, U, range, D, E, F are all not defined at that point.

In order to have a hope of optimizing your code, we need to be able to run your code.

You should also indicate which of the variables are always going to be the same for this purpose, and which you might change in the future.

Inside function1:

 spektras1 = mspd1; %initial file 1
 spektras2 = mspd2; %iniitial file 2
 spektras3 = mspd3; %initial file 3

None of those mspd* variables are defined.

http://www.mathworks.com/help/matlab/math/parameterizing-functions.html

[Domantas B]

This is my main code, I will upload all functions which are used.

clc
clear all

tic

range = 451:680;
unitconversion = (1e-6/.0001);

Xn = 2.4642/2.4830*100;
Yn = 1*100;
Zn = 1.8828/2.4830*100;

infile1 = 'A.txt';
mspd1 = load(infile1, '-ascii');

infile2 = 'B.txt';
mspd2 = load(infile2, '-ascii');

infile3 = 'C.txt';
mspd3 = load(infile3, '-ascii');

infile4 = 'D.csv';
mspd4 = load(infile4, '-ascii');

infile5 = 'E.dat';
mspd5 = load(infile5, '-ascii');

xconst = 0.3608;
yconst = 0.3635;

[ anew1, ctrwat, dcbat ] = function3(mspd4);

 b = 0;
 ci = 0;
 cj = 0;
 ck = 0;
 kr1 = 0;
 kg2 = 0;
 kb3 = 0;
 
 
for i = -20:1:20
    for j = -20:1:20
        for k = -20:1:20
            [spd,kr, kg, kb] = function1(i,j,k,mspd1,mspd2,mspd3);   
            deltabendras = function2( mspd4,mspd5,Xn,Yn,Zn,unitconversion,range,anew1,ctrwat,dcbat,spd );
            if(b<=deltabendras) 
                b = deltabendras;
                ci = i; cj = j; ck = k;
                kr1 = kr; kg2 = kg; kb3 = kb; 
            end
        end
    end
end

[Domantas B]

https://drive.matlab.com/sharing/b62356cd-55f5-482a-b7a6-8259fae1c70e

I put all my functions and files here, because mathworks limit uploads to 10 files.

Now the code should work.

Answer 1

[Walter Roberson]

syms kr kg kb
eqn1 = raudonaXnormuotas*kr + zaliaXnormuotas*kg + melynaXnormuotas*kb == xconst;
eqn2 = raudonaYnormuotas*kr + zaliaYnormuotas*kg + melynaYnormuotas*kb == yconst;
eqn3 = kr + kg + kb == 1;
eqns = [eqn1 eqn2 eqn3];

sol = solve(eqns, [kr, kg, kb]);
kr = double(sol.kr);
kg = double(sol.kg);
kb = double(sol.kb);

You do not need the symbolic toolbox at all for that code. You can use

A = [
    raudonaXnormuotas, zaliaXnormuotas, melynaXnormuotas;
    raudonaYnormuotas, zaliaYnormuotas, melynaYnormuotas;
                    1,               1,                1
    ]
b = [xconst; yconst; 1];
kvals = A\b;
kr = kvals(1);
kg = kvals(2);
kb = kvals(3);

[Walter Roberson]

rangespektras = spektras1iki471(92:1:321,:);
outspektras1 = spektras1iki471(kieknm+92:1:kieknm+321,:); %praleidzia pirmus 50nm
outspektras2 = spektras2iki471(kieknm2+92:1:kieknm2+321,:); %praleidzia pirmus 50nm
outspektras3 = spektras3iki471(kieknm3+92:1:kieknm3+321,:); %praleidzia pirmus 50nm
outspektras1nm = [outspektras1(:,1)-kieknm, outspektras1(:,2)];
outspektras2nm = [outspektras2(:,1)-kieknm2, outspektras2(:,2)];
outspektras3nm = [outspektras3(:,1)-kieknm3, outspektras3(:,2)];


XYZ2degred = gettristimulus2deg(interpolatespd([outspektras1nm(:,1) outspektras1nm(:,2).*unitconversion],rangered),rangered);
xyz2degred = getxyz(XYZ2degred);

raudonaX = XYZ2degred(1,1);
raudonaY = XYZ2degred(1,2);
raudonaZ = XYZ2degred(1,3);

raudXYZsuma = raudonaX+raudonaY+raudonaZ;

You are re-doing that calculation for every i, j, k combination, but it only changes with i.

You should create a routine that calculates those given i (and whatever appropriate inputs), and you should call that just inside your for i loop, getting the result, and passing it into function1.

XYZ2deggreen = gettristimulus2deg(interpolatespd([outspektras2nm(:,1) outspektras2nm(:,2).*unitconversion],rangered),rangered);
xyz2deggreen = getxyz(XYZ2deggreen); 

zaliaX = XYZ2deggreen(1,1);
zaliaY = XYZ2deggreen(1,2);
zaliaZ = XYZ2deggreen(1,3);

zaliaXYZsuma = zaliaX+zaliaY+zaliaZ;

That only changes with j (and not with i or k). You should create a routine that calulates those given j (and whatever appropriate inputs), and before your for i loop, you should calculate those for all possible j values and store them, and then just inside the for j loop you should pull out the appropriate one and pass it in to function1.

XYZ2degblue = gettristimulus2deg(interpolatespd([outspektras3nm(:,1) outspektras3nm(:,2).*unitconversion],rangered),rangered);
xyz2degblue = getxyz(XYZ2degblue);

melynaX = XYZ2degblue(1,1);
melynaY = XYZ2degblue(1,2);
melynaZ = XYZ2degblue(1,3);

melynaXYZsuma = melynaX+melynaY+melynaZ;

Changes only with k, create a routine, store the values, yada yada yada.

You should avoid re-calculating anything substantive. Anything like this that can be pre-calculated, should be pre-calculated.

[Domantas B]

Thanks a lot!

When I only changed that kr,kg,kb calculations to without symbolic toolbox, my code now is 20-30 times faster. I will try to avoid re-calculating to make better results.