نحوه حرکت روی خط در متلب

شما می دونید که معادله یک خط در فضا به شکل y=ax+b هستش و  نمایش دیگر و کاربردی اون هم به شکل (y-y1=m(x-x1  هستش، حالا در اینجا ما فقط یک خط در دو بعد داریم، با این حساب ب داشتن شیب یا همان m  و یک نقطه (x1,y1) ما قادر خواهیم بود که با گرفتن مقدار مختصه x مقدار تابع را به ازای آن مختصه بدست بیاوریم. یا برعکس با داشتن مقدار y شما می‌توانید مقدار x  نقطه را به ازای آن بدست بیاوید . در زیز کد متلب این عمل را می آوریم.

%% clear commands
delete(allchild(0));% delete all figures and graphical objects
clear;% Remove items from workspace, freeing up system memory
clc;% clear the Console


m=5;
x1=5, y1=10;

syms x;
y=m*(x-x1)+y1;
x2=input('give a x coordinate:\n ');
x=x2;
y2=eval(y);
plot([x1,x2],[y1,y2],'o--g');hold on;
 
x3=input('give another one x coordinate:\n ');
x=x3;
y3=eval(y);
plot([x2,x3],[y2,y3],'o--m');

نمونه خروجی این کد:

 

 

بخش اول کد پیاده سازی الگوریتم مورد استفاده برای رمزگذاری در bitcoin

همانطور که می‌دونید معادله عمومی منحنی بیضوی به شکل y² = x³ + ax + b می‌‍باشد که در ورژن بیت کوین مقدار ضریب a برابر صفر و مقدار ثابت b برابر 7 می‌باشد. نمای کلی منحنی بیضوی بیت کوین به شکل زیر می‌باشد:

 

از ویژگی‌های جالب این منحنی عبارتند از:
1-  اگر دو نقطه را در روی منحنی انتخاب کنید که خط واصل آن‌ها عمودی نباشد، خط واصل این دو نقطه، نقطه سومی را روی منحنی قطع خواهد کرد، که اگر این  نقطه را نسبت به محور xها متقارن کنیم، حاصل آن برابر با حاصل جمع دو نقطه اولیه خواهد بود...(point adding)

برای اضافه کردن دو نقطه به هم از روابط زیر استفاده می‌شود، که rx و ry حاصل اضافه کردن این دو نقطه می‌باشد.

c = (qy – py) / (qx – px)
rx = c2 – px – qx
ry = c (px – rx) – py

2- اگر یک نقطه ای را در روی این نمودار انتخاب کنیم، و چند بار ان نقطه را با خودش جمع کنیم، دوباره نقطه جدید روی این منحنی خواهد افتاد.(point doubling)

برای اضافه کردن یک نقطه به خودش از روابط زیر استفاده می‌شود، که rx و ry حاصل اضافه کردن این نقطه به خودش می‌باشد.

c = (3px2 + a) / 2py
rx = c2 – 2px
ry = c (px – rx) – py

 

در ادامه کد متلب پیاده‌سازی این عملیات را می‌توانید به شکل زیر مشاهده کنید:

clc;
clear;
close all;
% An elliptic curve is represented algebraically as an equation of the form
% y2 = x^3 + ax + b  ,For a = 0 and b = 7 (the version used by bitcoin), it looks like this y2 = x^3 + 7 

x=- nthroot(7,3):0.02: nthroot(7,3)*5; y = sqrt(((x.^3) +7));
 plot(x, real(y),'r- ', 'lineWidth',1);hold on;
 plot(x,-real(y),'r- ', 'lineWidth',1);hold on;



%% point adding
px=x(5) ;py=real(y(5)) ;qx=x(500) ; qy=real(y(500)) ;
c = (qy - py) /  (qx - px) ;
rx = (c^2)  - px - qx;
ry = c* (px-rx) -py;
 [X,Y]=polyxpoly( [px,qx],[py,qy],x,real(y)); 
plot([px,qx],[py,qy],'r*-')
plot( X,Y,'go-');hold on;
plot( rx,ry,'ko-');hold on;
text( rx,ry,'point adding');


%% point doubling
px=x(60) ;py=real(y(25)) ;
c = (3*(px^2) + 0) / (2*py);
rx = (c^2)  - (2*px);
ry = c*(px - rx) - py;

plot( px,py ,'cd-')
plot( rx,ry,'co-');
text( rx,ry,'point doubling');

نمونه خروجی:

این پست با پیاده‌سازی کامل بروزرسانی خواهد شد.....

مطالعه بیشتر...

کد تبدیل زوایه دو قطبی به یک قطبی در متلب

همانطور که می دونید متلب زوایه 360 درجه را به شکل دو بخش منفی و مثبت 180 درجه ای نشون می ده، این مسئله باعث بروز مشکل در پیاده سازی برخی الگوریتم‌هایی می شه که با زوایه کار می کنند، به عنوان مثال زمانی که شما تعدادی نقاط پخش شده در فضا دارید و بخواهید که کوتاهترین مسیر بین آنها را پیدا کنید، یک روش ساده مرتب سازی این نودها بر مبنای مختصات قطبی و خصوصا زوایه هستش، در چنین مواقعی نقاط قطبی دو بخشی پاسخگو نیست بر همین اساس من یک تابع ساده را برای حل این شکل نوشتم، که ورودی اون زوایه دو بخشی(0تا 180 و 0تا -180) و خروجی اون زوایه یک بخشی 0تا 360 درجه هستش.

کد متلب این تابع:

 

%convert 2 pole angle  to 1 pole degree or radian  
%angl=Angle2pol_To_1pol([5,95,-95,0],2)
%angl=Angle2pol_To_1pol( -95,2)
function angl=Angle2pol_To_1pol(angl,type)
if nargin==1
type=1;
end
switch type
    case 1% radian type
        angl=mod(angl,2*pi);
        miInd=angl<0;
        angl(miInd)=(2*pi)+angl(miInd);
        
    case 2
        angl=mod(angl,360);
        miInd=angl<0;
        angl(miInd)=(360)+angl(miInd);
end
end

نمونه استفاده از این تابع و خروجی حاصل در  تصویر زیر ابتدا مختصات دکارتی نکات درهم و برم به مختصات قطبی ان تبدیل شده و سپس با تبدیل زوایه دو بخشی به زاوایه یک بخشی 360 درجه و مرتب سازی آن بر اساس زاویه، شکل سمت چپ به شکل سمت راست تغییر پیدا کرده، این روش بسیار ساده می تونه برای پیدا کردن کوتاهترین مسیر هم مورد استفاده قرار بگیره.

 

کد یک ویلای نیمه کاره در متلب

 

%%%%%%%%%%%%%%%%%%%%%%%%%%-----------simulation area-----------%%%%%%%%%%%%%%%%%%%%%%
% devloped in Iran(ardebil)                                 %
% Copyright (C) 2017(1396) by moussa hasanzadeh.           %
%All rights reserved.                                      %
%programed in mathlab 2016                                %
%gmail:[e-mail]
%phone: [phone]


function vlia1
        clc
        clear;
        delete(allchild(0));
        sideN=4;
        t1=linspace(0,1,sideN+1);%vector between 0 and 1
        X=-sin(2*pi*t1+(pi/sideN));%node X coordinates with phaze shift
        Y=-cos(2*pi*t1+(pi/sideN));%node Y coordinates with phaze shift
        X=((((X-min(X))/(max(X)- min(X)))*2)-1)./5;%normalize  and schale X
        Y=((((Y-min(Y))/(max(Y)- min(Y)))*2)-1)./5;%normalize and schale Y
        Z= (1:length(Y))*0 ;

        for x=0:20:40
            for y=0:20:40
                for z=1:15
                    plot3(X+x,Y+y,Z+z,'k');hold on
                    plot3([X+x;X+x],[Y+y;Y+y],[Z+z*0;Z+z],'k');hold on
                end
            end
        end



        xVec=[0,20,40,40,20,0,0];
        yVec=[0,0,0,40,40,40,0];
        zVec1=yVec*0;
        zVec1([2,5])=zVec1([2,5])+3;
        [xVec,yVec]=dividerE(xVec,yVec,10);

        zVec=yVec*0;
        [~,zVec1]=dividerE(xVec,zVec1,10);



        hj1=fill3([xVec],[yVec],[zVec+5],'w');alpha(hj1 ,1)
        hj1=fill3([xVec],[yVec],[zVec+5+0.2],'k');alpha(hj1 ,0.2)

        hj3=fill3([xVec;xVec],[yVec;yVec],[zVec+5;zVec+5+0.2],'k');alpha(hj3 ,0.7)

        hj2=fill3([xVec],[yVec],[zVec+10],'w');alpha(hj2 ,1)
        hj2=fill3([xVec],[yVec],[zVec+10+0.2],'k');alpha(hj2 ,0.2)

        hj3=fill3([xVec;xVec],[yVec;yVec],[zVec+10;zVec+10+0.2],'k');alpha(hj3 ,0.7)

        hj3=fill3([xVec],[yVec],[zVec+15],'k');alpha(hj3 ,0.2)
        hj3=fill3([xVec;xVec],[yVec;yVec],[zVec+15;zVec+15+0.2],'k');alpha(hj3 ,0.7)

        xVecPentHous=[0,20,40,20,0];
        yVecPentHous=[0,0,0,0,0] ;
        zVecPentHous=[15,15,15,18,15];
        for mm=1:39
            % fill3(xVecPentHous,yVecPentHous+mm,zVecPentHous,'k--')
            km=fill3(xVecPentHous,yVecPentHous+mm,zVecPentHous,'g-');hold on
            alpha(km,0.2)
        end

        fill3([xVec],[yVec],[zVec+15],'r')



        zVec1=zVec1+15;

        fill3([xVec],[yVec],[zVec1],'r-','edgecolor','k')
        fill3([xVec],[yVec],[zVec1-0.2],'r-','edgecolor','k')

        hj1=fill3([xVec*2]-20,[yVec*2]-20,[zVec],'g');alpha(hj1 ,1)
        yardVect=[13:length(xVec),1:10]

        hj1=fill3([[xVec(yardVect)*2]-20;[xVec(yardVect)*2]-20],[[yVec(yardVect)*2]-20;[yVec(yardVect)*2]-20],[zVec(yardVect);zVec(yardVect)+2],'r');




        for difIndex=1:10
            zVec=zVec+0.2;
            hj1=plot3([xVec(yardVect)*2]-20,[yVec(yardVect)*2]-20,[zVec(yardVect)],'ks-','lineWidth',1);alpha(hj1 ,1);hold on

        end
        zVec(11:12)=zVec(11:12)+1;
        xVec(11)= xVec(10);
        xVec(12)= xVec(13);
       
        hj1=plot3([xVec*2]-20,[yVec*2]-20,[zVec],'k-','lineWidth',2);alpha(hj1 ,1);hold on
         zVec([10,13])=0;
        
         door=[10,11,12,13,10]
         matrixObj.x=[xVec(door)*2]-20;
         matrixObj.y=[yVec(door)*2]-20;
         matrixObj.z=[zVec(door)];
         [ newx, newy, newz]=rotateParameter(matrixObj,[0,0,1],90);
         
         hj1=fill3(newx-sqrt((xVec(10)-xVec(13)).^2)*1.3, newy-32, newz ,'w','lineWidth',2);alpha(hj1 ,1);hold on

        

        end


        function   [X,Y]=dividerE(X,Y,n)

        T1=0;
        T2=0;

        for i=1:length(X)-1
            temp1=0;
            temp1=linspace(X(i),X(i+1),n+2)  ;

            T1=[T1,temp1(1:end-1)];
        end
        T1=[T1,temp1(end)];
        X=T1(2:end);


        for i=1:length(Y)-1
            temp2=0;
            temp2=linspace(Y(i),Y(i+1),n+2)  ;

            T2=[T2,temp2(1:end-1)];
        end
        T2=[T2,temp2(end)];
        Y=T2(2:end);

        end






            function [ newx, newy, newz]=rotateParameter(matrixObj,azel,alpha)

            u = azel(:)/norm(azel);
            alph = alpha*pi/180;
            cosa = cos(alph);
            sina = sin(alph);
            vera = 1 - cosa;
            x = u(1);
            y = u(2);
            z = u(3);
            rot = [cosa+x^2*vera x*y*vera-z*sina x*z*vera+y*sina; ...
                x*y*vera+z*sina cosa+y^2*vera y*z*vera-x*sina; ...
                x*z*vera-y*sina y*z*vera+x*sina cosa+z^2*vera]';

            x = matrixObj.x;
            y = matrixObj.y;
            z = matrixObj.z;

            [m,n] = size(z);
            if numel(x) < m*n
                [x,y] = meshgrid(x,y);
            end

            [m,n] = size(x);
            newxyz = [x(:) , y(:) , z(:) ];
            newxyz = newxyz*rot;
            newx =   reshape(newxyz(:,1),m,n);
            newy =  reshape(newxyz(:,2),m,n);
            newz =  reshape(newxyz(:,3),m,n);

            end

این هم یک تصویر از نتیجه کار: